iOS のUU数を調べる
背景
クッキーが全て取得できない場合でも UU 数を知りたい。
iPhone などのiOSのブラウザでは第三者クッキー(別ドメインのクッキー)が焼けないようになっている。 これは、広告を運用する上では非常に大きな問題である。 なぜなら広告のドメインは広告が表示されるページのドメインではないから、第三者クッキーとなるからである。
今ある手がかりを元になんとか統計を使ってUU数を調べる方法はないのだろうか。
方法
早速始めよう。
用意するデータ
同じ UU がアクセスした頻度を imps とそのUU数を表したデータを準備する。 例えば、一回しかアクセスしていない UU は imps = 1 の UU数に、 2回アクセスした場合には imps = 2 の UU 数にカウントする。
例えば次のようになる。
imps | uu |
---|---|
1 | 141497411 |
2 | 1952850 |
3 | 929600 |
4 | 677190 |
5 | 479118 |
6 | 377811 |
7 | 301407 |
8 | 252147 |
9 | 212153 |
10 | 184203 |
11 | 159196 |
12 | 141028 |
13 | 123652 |
14 | 111885 |
15 | 101167 |
16 | 91964 |
17 | 84616 |
18 | 78076 |
19 | 71527 |
20 | 66685 |
... | ... |
上記のデータは、ログから次のSQLで取得するイメージである。
select imps, count(1) uu from ( select uuid, count(1) imps FROM log_data group by uuid) group by imps
面倒なことに、クッキーが焼けていない UU は全て imps = 1 に混じってカウントされいる。これはクッキーの仕様上仕方ないこと。 (クッキーが焼けた場合のみカウントしようとすると、2回目以降でしか確認できないため、imps = 1 は原理上カウントできないことになってしまうため)
この imps = 1 の中に焼けていないクッキーも含まれていて、この焼けていないクッキーのUU数と焼けているがimps = 1 のUU数を推定するのが、今回の課題となる。
べき乗則を利用する
であるが、今回のクッキーの話に適用する。
先の頻度 imps ごとの UU 数を、imps をで、そのUU数を で表すとする。すると、
の式がべき乗則に従った分布の場合成り立つと考えられる。( はなんからの値)
数式だけでは分かりにくいと思うので、先のデータの impsとUU数の関係を二次元のグラフで見てみよう。
縦軸と横軸は対数にしてあるが、頻度の対数とUU数の対数のプロットが直線上に並んでいる。この様にプロットが直線上に並ぶときべき乗則に従った分布であると考えられる。
この直線を imps = 1 まで伸ばしていけば、頻度 imps = 1 の UU 数が推定できる。それには、この直線(線形関数)を決めているパラメータ を求める必要がある。
この を求める方法は色々あるが、今回は最小二乗法で求めることにしたい。
線形関数の最小二乗法は解析的に求めることができるので特別な分析ソフトウエアなどは必要としない。 ExcelやSQL上で容易に計算することができる。
例えば、SQLで求めるなら次のような感じである。
SELECT (n*d-c*e)/(n*b-e*e) AS a, (b*c-d*e)/(n*b-e*e) AS b, sum_uu FROM (SELECT sum(ln(imps)*ln(imps)) AS b, sum(ln(uu)) AS c, sum(ln(imps)*ln(uu)) AS d, sum(ln(imps)) AS e, count(1) AS n, sum(uu) AS sum_uu FROM imps_uu)
この例の場合、
が得られる。
これを使うと
だから、
と、クッキーが焼けて頻度が imps = 1 のUU数は、4927084 UUと推定できる。
クッキーが焼けていない UU 数
残るはクッキーが焼けていない UU 数の推定である。クッキーが焼けていないデータは、元のデータの imps =1 のUU数から、先の推定で分かった UU を引いたものである。 この例では、
が、そのインプレッションとなる。これはUU数ではないことに注意する。
さて、このインプレッションがどのようなUU数になるかは、クッキーが焼けている既知の頻度とUU数の分布に従うとしよう。 全インプレッション数がすでに求まっているのであれば、クッキーが焼けている既知のインプレッション数は、全インプレッション数からクッキーが焼けていないインプレッション数を引けば良い。 この例だと、全インプレッション数は、361042790 であるので、
である。これを使って残りのインプレッションがどこに分布するか計算できる。 その分布にしたがった残りのインプレッションの各頻度のインプ数からその頻度で割ることでUU数が求まる。 それらの合計のUU数がクッキーが焼けていないUU数である。すなわち、次の式で表せる。
ここで、 は、それぞれ、 残りのインプレッション数、クッキーが焼けているUUのインプレッション数である。 この例だと、 である。 最終的に、クッキーが焼けていない UU 数は、8199399 となった(少数切り捨て)。
全UU数
クッキーが焼けているUU数+クッキーが焼けていない UU 数となる。
クッキーが焼けているUU数は、imps = 2 以上の UU数と、べき乗則によって推定した imps = 1 のUU数である。 この例の場合、
である。また、クッキーが焼けていない UU 数は8199399 だった。
よって、求める全UU数は、
である。
まとめ
クッキーが取得できない場合に、既に取得できているクッキーの頻度データにもとづき、UU数を推定する手法について述べた。特殊な計算は必要なくSQLやExcelだけで計算できることを示した。
この記事は、VOYAGE GROUP Advent Calendar 2017 - Qiita に参加しています。
参考資料
今回使ったデータ(Excel の計算も有り)
参考サイト
ある確率未満であることを確認するのに必要なサンプルサイズ
背景
ある確率未満であることを確認するのに、どれだけサンプルサイズが必要だろうか。 例えば、0.5% 以上のクリック率(CTR)の広告枠にしか広告を出したくないという話を考えよう。 一度その広告枠に広告を出せば、0.5% 未満であること自体は信頼上限の計算にて確認できるが、毎日確率が変動する場合、ある日急に0.5%以上になることも考えられるので、たとえ0.5%未満であっても毎日広告を出し続けなければならない。 ではこの場合最低何回広告を表示させれば、クリック率が0.5%未満であるか否かが計算できるのであろうか。
導出
信頼上限の計算から逆算すれば、信頼上限がある確率未満になるのに必要なサンプルサイズが求められそうである。
Wilson score interval にて信頼上限を求めるとした場合を考える。
N
を全体の試行回数(求めたいサンプルサイズ)、 X
を成功回数、Z
を正規分布の分位数、P
をこれ未満にしたい信頼上限とする。このとき、Wilson socre interval の信頼上限は次のようになる。
P = (X+Z*Z/2+Z*sqrt(X*(N-X)/N+Z*Z/4))/(N+Z*Z)
これを maxima で解くことを考える。 以下やってみた。
(%i5) solve((X+Z*Z/2+Z*sqrt(X*(N-X)/N+Z*Z/4))/(N+Z*Z)=P,[N]); 2 2 N Z - 4 X + 4 N X 2 2 Z sqrt(-------------------) - 2 P Z + Z + 2 X N (%o5) [N = -----------------------------------------------] 2 P
解こうとしたが、sqrt が外せないらしい。そこで危険ではあるが私が式を最小限変形して(過去に何度も間違えた経験あるので載せることにする)、解ける形にする。
左辺の分母でかけて、
左辺の左の2つの項を右辺へ移動して、二乗根を取り除ける体制に、
慎重にやったので問題ないだろう。 左辺の二乗
Z^2*(X*(N-X)/N+Z*Z/4)
から、右辺
P*(N+Z*Z)-(X+Z*Z/2)
を引いたら0という形で、maxima の solve()
に再びかける。
(%i13) display2d:false; (%o13) false (%i14) solve(%o6-%o7^2,[N]); (%o14) [N = -Z^2,N = -(Z*sqrt((P^2-2*P+1)*Z^2+(4-4*P)*X)+(P-1)*Z^2-2*X)/(2*P), N = (Z*sqrt((P^2-2*P+1)*Z^2+(4-4*P)*X)+(1-P)*Z^2+2*X)/(2*P)]
2つの解が出たが、サンプルサイズは正の数であり、後者だろう。よって以下が求める式である。
ところで、X
はどうやって決めればよいのだろうか。普通の解釈をすれば、X
は動的な値である。すなわち、サンプルの取得中に、X
は増えていくので、それに応じて N
を増やすような形になる。
ところが、多くの場合は動的に取得できることは少ないので、サンプル数は増えるかもしれないが前回の値を直接使えばいいだろう。これは、信頼上限をこえた場合は一部のサンプルではなく、全取得するという前提である。
シミュレーション例
python で簡単なものを書いてみた。クリック率 0.5% 未満には広告を出したくないが、0.5%未満であるかどうかの探索はしたいという例である。
import math def sample_size(p, x, z=1.96): return (z*math.sqrt((p**2-2*p+1)*z**2+(4-4*p)*x)+(1-p)*z**2+2*x)/(2*p) def wilson_ucb(n, x, z=1.96): return (x+z*z/2+z*math.sqrt(x*(n-x)/n+z*z/4))/(n+z*z)
0クリックの枠
# 1 クリックもなかった広告枠に必要なサンプルサイズは sample_size(0.005, 0)
764.4784
# そのサンプルサイズで0クリックだった場合のCTRの信頼上限は wilson_ucb(765, 0)
0.004996607883860603
信頼上限が 0.5% を超えることはなかった。再び配信数(インプレッション)は制限される。
CTR が0.5% よりずっと低い枠
# 1万インプレッションで、5 クリックあった場合、CTR=0.05%だが、この枠のサンプルサイズは? sample_size(0.005, 5)
2336.485083416135
# 2337 インプレッション取得したところ、1クリックあった。CTRの信頼上限は、 wilson_ucb(2337, 1)
0.002068315942437113
信頼上限が、0.5% を超えることはなかった。再び配信数(インプレッション)は制限される。
# 1クリックだったのでサンプルサイズは? sample_size(0.005, 1)
1129.0503977791236
# 0 クリックだったのでサンプルサイズは? sample_size(0.005, 0)
764.4784
# 以上、クリックがなければ上記か、クリックがあれば、2つ上を繰り返す。
CTR が 0.5% 付近の場合
# 1000インプレッションで、5 クリックあった場合。 sample_size(0.005, 5)
2336.485083416135
wilson_ucb(1000, 5)
0.011242912497455987
# 信頼上限は1.1%。この時点で、全取得と思われるが、あえてそれをせずシミュレーション # 2337 インプレッションで、11 クリックあった場合。 sample_size(0.005, 11)
3934.259645857473
徐々に増えていくことがわかる。いずれ、その広告枠のインプレッション数の最大値に到達することだろう。
結論
もはや適切なサンプルサイズを調べるために、当選確実の計算方法や、統計のページを検索する必要はない。
N = (Z*sqrt((P^2-2*P+1)*Z^2+(4-4*P)*X)+(1-P)*Z^2+2*X)/(2*P)
が、必要なサンプルサイズだ。ここで、N
を全体の試行回数(求めたいサンプルサイズ)、 X
を成功回数、Z
を正規分布の分位数、P
をこれ未満にしたい信頼上限とする。
scikit learn の 確率推定のクラスの順序はソートされている。※ただしバイナリクラスの重みは…
背景
liblinear libsvm の重みや、確率推定の配列は、訓練データのクラスラベルの出現順で、決まるという恐ろしい仕様。 出現順というのは、1行目の事例のほうが2行目の事例より前、という出現順。 もちろん学習後、その順序が参照できるように、モデルの出力ファイルにはその出力クラスラベルの順序が書かれてるのだが、 それを読み取るのが面倒だったりするので、学習前に訓練データをクラスラベルでソートしたり tips があったりする[要出典]。
では liblinear libsvm を使っている、scikit learn はどうか?
調査
次のドキュメントを読むと、
sklearn.linear_model.LogisticRegression — scikit-learn 0.19.1 documentation
self.classes_
にあるらしい。
そして、 - numpy の unique() https://github.com/scikit-learn/scikit-learn/blob/a2ebb8cfd2d126ad8e6fb36e0bdadba7de8fcd9f/sklearn/preprocessing/label.py#L96
- 組み込みの sorted() https://github.com/scikit-learn/scikit-learn/blob/a2ebb8cfd2d126ad8e6fb36e0bdadba7de8fcd9f/sklearn/utils/multiclass.py#L105
もっと調べた
なんと、バイナリクラスの場合、逆転させていた。
liblinear 内の学習では、0 以上かどうかで判定している。
では、この prob_col->y
はいつどこで決められているのだろう。
prob->y
にて。
Y
→ y_ind
y
これが入力。
よって、次の fit_transform
がポイントになるが、
ためしに、[2, 2, 0, 0, 2]
を unique
で呼び出してみた。
>>> import numpy as np >>> np.unique([2, 2, 0, 0, 2], return_inverse=True) (array([0, 2]), array([1, 1, 0, 0, 1]))
すなわち、二値(バイナリ)の場合は後者のラベルが、prob_col->y
で0以上となることから、後者の重みを学習していることが分かる。
結論
sklearn の係数重みの出力順は、ラベルのソート順(昇順)になっている。
ただし、バイナリクラスの [0, 1] の重みを取ってくる場合は、coef_[0]
や、intercept_[0]
を取ってくるとラベルが1の方の重みである。
Logistic Regression における Negative Down Sampling
negative down sampling をしたときの logistic regression の最適な方法・設定を調べる。
結論としては、negative down sampling をする前に、 を平均から求める( )。それを対数化すれば切片の重みとなる( )。negative down sampling 後のデータは class_weight='balanced', fit_intercept=False
で学習を行う。
背景
logistic regression では サンプルの割合を変化させると学習はうまくいかない。 http://quinonero.net/Publications/predicting-clicks-facebook.pdf では、negative down sampling を行っているが、Re-Calibration という計算において、down sampling により曲げられた予測から再補正を行うことにより、もとの予測に戻している。 ところが、log-loss を損失関数にして考えた場合で、かつそれを元に学習パラメータの最適化を行った場合、この補正では満足な結果が得られない。理由としては、down sampling の前と後では、log-loss の値が変わり、そして、その極大点(損失関数的には neg_log_loss なので極小点)も異なるところにあると考えられるからである。 検証では次のサイトの例を参考にした。 https://stats.stackexchange.com/questions/67903/does-down-sampling-change-logistic-regression-coefficients
ここで確認したいこと
- intercept への罰則は適切ではないので、
intercept_scaling
を大きく設定 - negative down sampling することで neg_log_loss の極値となるパラメータに違いはあるか
class_weight
を使った補正の有効性- balanced の検証
intercept_scaling
を大きく設定したほうが良いか実験用データ生成
FAKE = 10 # フェイク変数の数
import pandas as pd import random def sample_data(size=10000, fake=FAKE): """ サンプルデータ生成。 x==1 のとき 10%, x==0 のとき 1%。 """ y = [] x = [] z = [[] for i in range(fake)] # fake r = [] for i in range(size): xx = 1 if random.random() < 0.05 else 0 x.append(xx) if xx == 1: yy = 1 if random.random() < 0.1 else 0 else: yy = 1 if random.random() < 0.01 else 0 y.append(yy) for j in range(fake): zz = 1 if random.random() > 0.1 else 0 z[j].append(zz) rr = random.random() if yy == 0 else 0 r.append(rr) data_dict = {'z'+str(i): z[i] for i in range(fake)} data_dict['y'] = y data_dict['x'] = x data_dict['r'] = r df = pd.DataFrame.from_dict( data_dict ) return df df_test = sample_data() df_test[df_test['r']< 0.005][:20]
r | x | y | z0 | z1 | z2 | z3 | z4 | z5 | z6 | z7 | z8 | z9 | |
---|---|---|---|---|---|---|---|---|---|---|---|---|---|
83 | 0.000000 | 0 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 |
132 | 0.000000 | 0 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 0 | 1 | 1 | 1 |
147 | 0.000000 | 0 | 1 | 1 | 1 | 1 | 1 | 1 | 0 | 1 | 1 | 1 | 1 |
324 | 0.005000 | 0 | 0 | 1 | 1 | 1 | 1 | 1 | 1 | 0 | 1 | 1 | 1 |
329 | 0.000000 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 |
363 | 0.001359 | 0 | 0 | 1 | 1 | 1 | 0 | 1 | 1 | 0 | 1 | 1 | 1 |
419 | 0.001078 | 0 | 0 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 |
494 | 0.000486 | 0 | 0 | 1 | 1 | 1 | 1 | 1 | 0 | 1 | 1 | 1 | 1 |
647 | 0.000000 | 1 | 1 | 0 | 1 | 1 | 0 | 1 | 1 | 1 | 1 | 1 | 1 |
698 | 0.003799 | 0 | 0 | 0 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 |
720 | 0.000000 | 0 | 1 | 1 | 1 | 1 | 0 | 1 | 1 | 1 | 1 | 1 | 1 |
737 | 0.000000 | 0 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 0 | 1 | 1 | 1 |
793 | 0.000000 | 0 | 1 | 1 | 0 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 |
797 | 0.000000 | 0 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 |
809 | 0.000000 | 0 | 1 | 1 | 1 | 1 | 1 | 0 | 1 | 0 | 1 | 1 | 1 |
925 | 0.000000 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 |
1035 | 0.000000 | 0 | 1 | 1 | 1 | 1 | 1 | 1 | 0 | 1 | 1 | 1 | 1 |
1093 | 0.000000 | 0 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 |
1110 | 0.000000 | 0 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 |
1145 | 0.000000 | 0 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 1 | 0 | 1 |
df_list = [sample_data() for i in range(10)] df = df_list[0]
ダウンサンプリング率の変更による重みの確認
from sklearn.linear_model import LogisticRegression from IPython.core.display import display, Markdown def check_coef_by_down_sampling(df_list, intercept_scale=1): for ratio in [1, 2, 5, 10, 20, 50, 100, 200]: display(Markdown("### ratio: %s" % str(ratio))) for i in range(0, 10): df = df_list[i] df_d = df[df['r'] < 1.0/ratio] logreg = LogisticRegression(penalty='l1', C=1, intercept_scaling=intercept_scale) column_list = ['x'] + ['z'+str(i) for i in range(FAKE)] logreg.fit(df_d[column_list], pd.Series(df_d['y'].values.flatten())) w_x = logreg.coef_.tolist()[0][0] w_z = sum([w*w for w in logreg.coef_.tolist()[0][1:]]) # 二乗して合計(ここはフェイクなので小さい方がいい) print([round(x,5) for x in [logreg.intercept_[0], w_x, w_z] ])
check_coef_by_down_sampling(df_list)
ratio: 1
[-4.91204, 2.74654, 0.17877]
[-4.0515, 2.3485, 0.83755]
[-5.42401, 2.51179, 0.4848]
[-4.03337, 2.06947, 0.44515]
[-4.69719, 2.57411, 0.59647]
[-3.76437, 2.35255, 0.35005]
[-4.39169, 2.05972, 0.61103]
[-4.11071, 2.53458, 0.17263]
[-2.99168, 2.46254, 0.63564]
[-3.99944, 2.30797, 0.55562]
ratio: 2
[-4.2107, 2.70525, 0.23816]
[-3.38297, 2.31001, 0.92924]
[-4.77464, 2.50352, 0.55668]
[-3.37317, 2.03352, 0.4837]
[-4.03382, 2.53982, 0.58345]
[-2.97067, 2.31873, 0.38171]
[-3.41311, 2.12467, 0.67182]
[-3.44971, 2.49872, 0.19891]
[-2.39652, 2.40039, 0.5989]
[-3.3235, 2.31424, 0.65252]
ratio: 5
[-3.35169, 2.72117, 0.23848]
[-2.27844, 2.41117, 0.77865]
[-3.52712, 2.46208, 0.72024]
[-2.23162, 2.27092, 0.56094]
[-3.22692, 2.54749, 0.57494]
[-1.91845, 2.26757, 0.57179]
[-2.54464, 2.15011, 0.68264]
[-2.45929, 2.51433, 0.26709]
[-1.19256, 2.56356, 0.80125]
[-2.54425, 2.20229, 0.74773]
ratio: 10
[-3.25613, 2.75951, 0.52]
[-1.42508, 2.38864, 0.48791]
[-2.41798, 2.45274, 0.90139]
[-1.60434, 2.24338, 0.59875]
[-1.95903, 2.56728, 0.65356]
[-1.02283, 2.17213, 0.57539]
[-1.56087, 2.2855, 0.89761]
[-1.4418, 2.45697, 0.47024]
[-0.72239, 2.52637, 0.75559]
[-1.6296, 2.47578, 0.9896]
ratio: 20
[-2.40845, 2.86692, 0.62559]
[-1.03619, 2.22112, 0.40693]
[-1.43415, 2.61045, 0.5262]
[-0.45256, 2.03482, 0.52997]
[-1.19788, 2.51803, 0.47951]
[0.0, 2.08584, 0.65126]
[-1.23638, 2.23326, 1.23025]
[-0.46025, 2.78909, 0.70659]
[-0.45358, 2.60837, 0.5026]
[-1.25244, 2.28267, 0.91885]
ratio: 50
[-0.66989, 2.5969, 0.87268]
[0.0, 2.07957, 0.88407]
[0.0, 2.63922, 0.5973]
[0.0, 1.92209, 0.67016]
[-0.52315, 2.67275, 0.84047]
[0.0, 1.69736, 0.38221]
[0.0, 2.49584, 1.50871]
[0.0, 2.58153, 0.85577]
[0.0, 2.57994, 0.81271]
[-0.88678, 2.0582, 0.9175]
ratio: 100
[0.0, 3.11442, 0.84759]
[0.0, 1.93646, 1.44019]
[0.0, 2.61444, 1.39052]
[0.0, 1.8008, 0.48262]
[-0.67109, 3.00093, 1.18918]
[0.0, 1.37616, 0.1724]
[0.0, 2.52167, 1.31528]
[0.0, 2.44528, 0.80262]
[0.0, 2.23357, 0.57737]
[0.0, 2.59899, 0.89581]
ratio: 200
[0.0, 2.86792, 1.13424]
[0.0, 2.1795, 0.93253]
[0.2598, 2.56551, 1.55929]
[0.05154, 1.95667, 0.3683]
[0.0, 2.43416, 1.52139]
[0.0, 2.09143, 0.38268]
[0.0, 2.10479, 1.78338]
[0.18917, 2.81943, 0.89836]
[0.0, 3.27659, 1.14055]
[0.0, 2.41917, 0.91509]
一つ目は、切片、2つめは y
の生成に寄与している x
、それ以降は全く関係ない z*
の回帰係数の二乗の合計である。なお、二番目の x
の重みは、1%から10%へ10倍なので、10の自然対数である、2.302585 が理想的である。
ダウンサンプリング率によって
x
の回帰係数の重みはほとんど変化しないが、切片は大きく変わる。- ダウンサンプリング率が大きくなると、ノイズの
z*
の重みを誤って学習しやすい ということがわかる。
check_coef_by_down_sampling(df_list, intercept_scale=10000)
ratio: 1
[-5.45556, 2.75912, 0.28262]
[-4.6428, 2.35239, 0.92031]
[-5.78945, 2.52124, 0.64852]
[-4.49992, 2.07685, 0.49226]
[-5.20276, 2.58633, 0.69073]
[-4.34302, 2.36255, 0.28314]
[-5.01426, 2.0722, 0.7266]
[-4.6828, 2.54017, 0.13943]
[-3.47544, 2.47044, 0.48384]
[-4.54506, 2.31714, 0.56532]
ratio: 2
[-4.47498, 2.71107, 0.27932]
[-4.04952, 2.31372, 1.04347]
[-5.06899, 2.50929, 0.67163]
[-3.87872, 2.04466, 0.54536]
[-4.33338, 2.54452, 0.60214]
[-3.54863, 2.32828, 0.29453]
[-4.04781, 2.13091, 0.7376]
[-4.02237, 2.50532, 0.17511]
[-2.98376, 2.41174, 0.45404]
[-3.84655, 2.32399, 0.66282]
ratio: 5
[-3.97839, 2.73628, 0.37851]
[-2.93231, 2.40569, 0.85465]
[-4.28537, 2.47596, 0.97301]
[-2.35712, 2.27262, 0.55879]
[-4.03678, 2.56098, 0.76631]
[-2.53755, 2.27179, 0.45956]
[-3.13683, 2.15767, 0.73814]
[-2.9908, 2.51932, 0.22433]
[-1.8004, 2.56204, 0.57123]
[-3.09578, 2.21253, 0.79138]
ratio: 10
[-3.72013, 2.77988, 0.7003]
[-2.08104, 2.37616, 0.48348]
[-3.22456, 2.4633, 1.11463]
[-2.22255, 2.2543, 0.64177]
[-2.92231, 2.58176, 0.72124]
[-1.64388, 2.18194, 0.4027]
[-2.34855, 2.29581, 0.92035]
[-2.15464, 2.45468, 0.38139]
[-1.33011, 2.52442, 0.55185]
[-2.25828, 2.48406, 1.00449]
ratio: 20
[-3.00527, 2.89685, 0.8809]
[-1.62548, 2.2128, 0.43464]
[-2.2563, 2.62659, 0.68069]
[-0.98092, 2.04177, 0.4648]
[-1.47433, 2.52369, 0.45993]
[-0.57554, 2.09909, 0.42959]
[-1.85798, 2.26083, 1.32768]
[-1.37675, 2.77282, 0.49475]
[-1.26171, 2.60299, 0.35655]
[-1.74328, 2.28939, 0.97978]
ratio: 50
[-1.84861, 2.62367, 1.05701]
[-0.51107, 2.07379, 0.84463]
[-1.10631, 2.63535, 0.60667]
[-0.1167, 1.92674, 0.56876]
[-1.42434, 2.7109, 0.87041]
[-0.18693, 1.70358, 0.35005]
[-0.7797, 2.51998, 1.44281]
[0.09158, 2.58324, 0.88211]
[-0.46089, 2.57796, 0.76345]
[-1.71747, 2.1023, 1.21375]
ratio: 100
[-0.59101, 3.13524, 0.90419]
[-0.1334, 1.9426, 0.97734]
[0.07708, 2.61341, 1.3952]
[0.67081, 1.79631, 0.55785]
[-2.34283, 3.06153, 2.02438]
[-0.277, 1.39244, 0.18829]
[0.65186, 2.51285, 1.41825]
[0.3768, 2.44727, 0.88092]
[0.29966, 2.23463, 0.60469]
[-1.20692, 2.6587, 1.16563]
ratio: 200
[-0.40691, 2.90292, 1.28256]
[0.53202, 2.184, 0.87883]
[1.80984, 2.59295, 1.85358]
[0.95417, 1.92716, 0.19392]
[-0.92467, 2.48612, 1.94237]
[-0.38448, 2.11439, 0.50066]
[1.05276, 2.08541, 1.95901]
[2.15207, 2.8209, 1.84377]
[1.16178, 3.27734, 1.26292]
[-0.13707, 2.42296, 0.95436]
intercept_scaling を大きく設定する(後半の設定)
切片の重みが無理なく付くようになり、x
への重みも、2.302585 に近い値を安定して獲得できるようになった。 negative down sampling で意図的に変更している以上、切片の重みに罰則を与える意味は全く無い。
なお、切片の重みは、与えた intercept_scaling
ですでに乗算してあるらしく、乗算する必要はない。ソースはここ
https://github.com/scikit-learn/scikit-learn/blob/45dc891c96eebdb3b81bf14c2737d8f6540fabfe/sklearn/svm/base.py#L902
intercept_scaling
を大きく設定したほうが良い
- negative down sampling をしないほうがいいし、率は小さい方がいい。(それでもメモリの問題などの理由がある)
intercept_scaling
を大きく設定したほうが、無駄に切片の重みを学習せずにすむnegative down sampling することで neg_log_loss の極値となるパラメータに違いはあるか
探索用 cross validation
from sklearn.model_selection import GridSearchCV def grid_search_cv(X, y, class_weight=None, fit_intercept=True): test_parameters = [ {'C': [0.01, 0.1, 1, 10, 100] } ] clf = GridSearchCV( LogisticRegression(penalty='l1', intercept_scaling=10000, class_weight=class_weight, fit_intercept=fit_intercept), test_parameters, cv=20, scoring=['neg_log_loss', 'accuracy'], n_jobs=-1, refit='neg_log_loss' ) clf.fit(X, y) return clf
for ratio in [1, 10, 100, 200]: display(Markdown("### ratio: %s" % str(ratio))) df_d = df[df['r'] < 1.0/ratio] column_list = ['x'] + ['z'+str(i) for i in range(FAKE)] grid_search_result = grid_search_cv(df_d[column_list], pd.Series(df_d['y'].values.flatten())) df_result = pd.DataFrame(grid_search_result.cv_results_) display(df_result[['params', 'mean_test_neg_log_loss', 'std_test_neg_log_loss']])
ratio: 1
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.077886 | 0.004059 |
1 | {u'C': 0.1} | -0.068239 | 0.006209 |
2 | {u'C': 1} | -0.068588 | 0.007413 |
3 | {u'C': 10} | -0.068781 | 0.007660 |
4 | {u'C': 100} | -0.069028 | 0.007834 |
ratio: 10
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.390054 | 0.012941 |
1 | {u'C': 0.1} | -0.326330 | 0.028839 |
2 | {u'C': 1} | -0.326252 | 0.038349 |
3 | {u'C': 10} | -0.329250 | 0.041436 |
4 | {u'C': 100} | -0.329717 | 0.041785 |
ratio: 100
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.661588 | 0.005336 |
1 | {u'C': 0.1} | -0.576321 | 0.037786 |
2 | {u'C': 1} | -0.554997 | 0.121124 |
3 | {u'C': 10} | -0.571344 | 0.137919 |
4 | {u'C': 100} | -0.587997 | 0.181002 |
ratio: 200
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.554096 | 0.031750 |
1 | {u'C': 0.1} | -0.516677 | 0.032798 |
2 | {u'C': 1} | -0.473868 | 0.112612 |
3 | {u'C': 10} | -0.530726 | 0.213819 |
4 | {u'C': 100} | -0.554108 | 0.242809 |
negative down sampling することで neg_log_loss の極値となるパラメータに違いは確認できず
- neg_log_loss の値自体に違いがあった。(事例数が異なるから当たり前か)
- neg_log_loss の極値となる学習パラメータに多少の違いはあったが、誤差の範囲で、今回の検証では分からなかった。
class_weight
を使った補正の有効性class_weight
を使うことで失われた事例分を残った事例に重みに重ねることができる。これで学習時のバランスが補正されないかというアイデア。
import numpy as np def make_class_weight(orig, down, log): """クラスの重みを計算する。down sampling によって減った事例数を戻すように補正する""" neg_w = float(orig) / down if log: neg_w = np.log(neg_w + 1) return {0: neg_w, 1: 1} def check_learn_param(class_weight=False, log=False): for ratio in [1, 2, 5, 10, 20, 50, 100, 200]: display(Markdown("### ratio: %s" % str(ratio))) df_d = df[df['r'] < 1.0/ratio] cw = make_class_weight(len(df[df['y']==0]), len(df_d[df_d['y']==0]), log) if class_weight else None column_list = ['x'] + ['z'+str(i) for i in range(FAKE)] grid_search_result = grid_search_cv(df_d[column_list], pd.Series(df_d['y'].values.flatten()), class_weight=cw) df_result = pd.DataFrame(grid_search_result.cv_results_) display(df_result[['params', 'mean_test_neg_log_loss', 'std_test_neg_log_loss']]) check_learn_param(class_weight=True)
ratio: 1
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.077886 | 0.004059 |
1 | {u'C': 0.1} | -0.068239 | 0.006205 |
2 | {u'C': 1} | -0.068631 | 0.007474 |
3 | {u'C': 10} | -0.068918 | 0.007800 |
4 | {u'C': 100} | -0.069015 | 0.007842 |
ratio: 2
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.141263 | 0.008058 |
1 | {u'C': 0.1} | -0.122711 | 0.011923 |
2 | {u'C': 1} | -0.123334 | 0.014466 |
3 | {u'C': 10} | -0.123813 | 0.015084 |
4 | {u'C': 100} | -0.124233 | 0.015270 |
ratio: 5
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.310091 | 0.018207 |
1 | {u'C': 0.1} | -0.266433 | 0.026865 |
2 | {u'C': 1} | -0.267755 | 0.032059 |
3 | {u'C': 10} | -0.268953 | 0.034197 |
4 | {u'C': 100} | -0.269233 | 0.033790 |
ratio: 10
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.567181 | 0.029106 |
1 | {u'C': 0.1} | -0.485387 | 0.049120 |
2 | {u'C': 1} | -0.485032 | 0.061664 |
3 | {u'C': 10} | -0.486115 | 0.064668 |
4 | {u'C': 100} | -0.486117 | 0.063645 |
ratio: 20
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.991310 | 0.041022 |
1 | {u'C': 0.1} | -0.833435 | 0.086217 |
2 | {u'C': 1} | -0.830881 | 0.113298 |
3 | {u'C': 10} | -0.830402 | 0.123086 |
4 | {u'C': 100} | -0.833500 | 0.125014 |
ratio: 50
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -1.803340 | 0.068163 |
1 | {u'C': 0.1} | -1.551492 | 0.148318 |
2 | {u'C': 1} | -1.519491 | 0.213549 |
3 | {u'C': 10} | -1.532629 | 0.227115 |
4 | {u'C': 100} | -1.526625 | 0.225905 |
ratio: 100
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -2.630216 | 0.045231 |
1 | {u'C': 0.1} | -1.994720 | 0.314057 |
2 | {u'C': 1} | -1.969576 | 0.448479 |
3 | {u'C': 10} | -1.995761 | 0.481577 |
4 | {u'C': 100} | -1.961631 | 0.429366 |
ratio: 200
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -3.184818 | 0.122700 |
1 | {u'C': 0.1} | -2.228196 | 0.525339 |
2 | {u'C': 1} | -2.231107 | 0.676387 |
3 | {u'C': 10} | -2.318523 | 0.762673 |
4 | {u'C': 100} | -2.345855 | 0.750384 |
neg_log_loss
は予想に反して同じ値にならなかった。ただ、極値は保存されているようである。
係数の重みの確認
def check_coefs_for_each_ratio(df_list, intercept_scale=1, class_weight=False, log=False): for ratio in [1, 2, 5, 10, 20, 50, 100, 200]: display(Markdown("### ratio: %s" % str(ratio))) for i in range(0, 10): df = df_list[i] df_d = df[df['r'] < 1.0/ratio] cw = make_class_weight(len(df[df['y']==0]), len(df_d[df_d['y']==0]), log) if class_weight else None logreg = LogisticRegression(penalty='l1', C=1, intercept_scaling=intercept_scale, class_weight=cw) column_list = ['x'] + ['z'+str(i) for i in range(FAKE)] logreg.fit(df_d[column_list], pd.Series(df_d['y'].values.flatten())) w_x = logreg.coef_.tolist()[0][0] w_z = sum([w*w for w in logreg.coef_.tolist()[0][1:]]) # 二乗して合計(ここはフェイクなので小さい方がいい) print([round(x,5) for x in [logreg.intercept_[0], w_x, w_z] ])
check_coefs_for_each_ratio(df_list, intercept_scale=10000, class_weight=True)
ratio: 1
[-5.456, 2.75913, 0.28275]
[-4.68984, 2.35276, 0.93339]
[-5.93321, 2.52496, 0.73554]
[-4.49972, 2.07685, 0.49214]
[-5.07817, 2.58326, 0.65842]
[-4.18139, 2.36003, 0.29171]
[-5.01522, 2.07222, 0.72691]
[-4.63656, 2.53985, 0.13835]
[-2.93413, 2.46116, 0.66034]
[-4.58423, 2.31785, 0.56987]
ratio: 2
[-4.91084, 2.70796, 0.24419]
[-4.73899, 2.3222, 1.09332]
[-5.49024, 2.50968, 0.57039]
[-4.57879, 2.0469, 0.54657]
[-5.37927, 2.55492, 0.71975]
[-4.18683, 2.33364, 0.30944]
[-4.72929, 2.13681, 0.76036]
[-4.78652, 2.50797, 0.18094]
[-3.74175, 2.41311, 0.45017]
[-4.47882, 2.32848, 0.68062]
ratio: 5
[-5.52737, 2.75556, 0.44584]
[-4.32948, 2.43492, 0.95364]
[-5.36533, 2.48112, 0.79776]
[-4.49862, 2.30294, 0.67979]
[-4.66343, 2.56162, 0.58413]
[-3.80696, 2.3091, 0.61495]
[-4.68544, 2.18322, 0.82486]
[-4.55734, 2.52214, 0.21192]
[-3.49043, 2.57608, 0.64307]
[-4.55503, 2.22938, 0.85661]
ratio: 10
[-6.32183, 2.84784, 1.00543]
[-4.18039, 2.42215, 0.64773]
[-5.70742, 2.54972, 1.53197]
[-4.58312, 2.32635, 0.85909]
[-5.08476, 2.60348, 0.69908]
[-3.32753, 2.26653, 0.80517]
[-4.30632, 2.3637, 1.1924]
[-4.57177, 2.47107, 0.39637]
[-3.74324, 2.5844, 0.79115]
[-4.25131, 2.58967, 1.42814]
ratio: 20
[-3.73082, 2.97579, 1.38839]
[-4.48063, 2.30836, 0.82273]
[-5.26905, 2.81853, 1.29472]
[-4.26113, 2.13226, 0.80125]
[-5.23599, 2.62421, 0.66548]
[-2.60202, 2.29986, 1.61421]
[-4.92805, 2.39065, 2.24484]
[-4.5201, 2.85804, 0.65841]
[-4.47336, 2.67563, 0.54306]
[-4.70123, 2.38434, 1.21893]
ratio: 50
[-6.1262, 3.10075, 3.12108]
[-4.1213, 2.26967, 1.92844]
[-6.14338, 2.9595, 1.32346]
[-4.26593, 2.21141, 1.51868]
[-3.33964, 3.0612, 1.31777]
[-3.13537, 1.96733, 1.71331]
[-3.31786, 3.01146, 3.9546]
[-3.43166, 3.07311, 4.31343]
[-4.97486, 2.72593, 0.82405]
[-5.88227, 2.49085, 3.02388]
ratio: 100
[-8.0263, 5.20309, 14.74563]
[-3.96654, 2.50111, 5.92777]
[-3.77454, 3.32975, 4.06035]
[-2.66965, 2.29368, 2.19824]
[-9.43126, 5.73963, 11.59807]
[-4.36287, 1.76316, 1.47788]
[-1.57113, 3.0898, 4.08203]
[-4.62372, 3.20552, 5.05104]
[-4.9451, 2.44147, 0.9702]
[-8.34704, 4.66164, 11.21491]
ratio: 200
[-9.04834, 6.70274, 26.74032]
[-3.0757, 2.93908, 6.78269]
[4.31901, 3.79058, 15.57901]
[-3.1984, 2.82233, 2.86396]
[-9.31711, 5.73608, 16.86566]
[-9.31487, 5.56711, 19.90607]
[-0.56633, 2.71122, 15.52765]
[-6.83777, 6.20401, 20.82009]
[-1.64302, 8.69927, 9.08549]
[-7.99419, 5.54549, 9.48558]
切片も出ているが、x
への重みのばらつき(2.302585が理想)や、z
の重みへのばらつき(0が理想)が大きい。
重みの対数化
安直ではあるが、10倍するのではなく、log(10)倍することで、過学習を緩和したい
check_coefs_for_each_ratio(df_list, intercept_scale=10000, class_weight=True, log=True)
ratio: 1
[-5.06334, 2.75755, 0.2739]
[-4.285, 2.34939, 0.90332]
[-5.54231, 2.52186, 0.71662]
[-4.11011, 2.07534, 0.48425]
[-4.79435, 2.58338, 0.6687]
[-3.54648, 2.35229, 0.31755]
[-4.61396, 2.07039, 0.71198]
[-4.25854, 2.53891, 0.1378]
[-3.00963, 2.46912, 0.50412]
[-4.22915, 2.31719, 0.56841]
ratio: 2
[-4.9026, 2.71966, 0.36531]
[-4.14688, 2.31516, 1.05298]
[-5.48156, 2.51683, 0.83853]
[-3.98731, 2.04515, 0.54704]
[-4.81269, 2.55341, 0.71851]
[-3.62939, 2.32907, 0.29694]
[-4.14463, 2.13188, 0.74179]
[-3.78081, 2.50169, 0.1792]
[-3.0896, 2.41201, 0.45301]
[-3.6255, 2.31867, 0.64602]
ratio: 5
[-4.55296, 2.74552, 0.41394]
[-3.24339, 2.42041, 0.86243]
[-4.98751, 2.48683, 1.04366]
[-3.44419, 2.29147, 0.63844]
[-4.27659, 2.56317, 0.66561]
[-2.95918, 2.28709, 0.5224]
[-3.69191, 2.16886, 0.7741]
[-3.6336, 2.52211, 0.22204]
[-2.41829, 2.56805, 0.60344]
[-3.38312, 2.21611, 0.78615]
ratio: 10
[-4.97196, 2.82519, 0.98587]
[-2.9001, 2.39795, 0.56099]
[-3.63711, 2.48865, 1.09703]
[-3.22029, 2.29187, 0.7576]
[-3.71282, 2.59307, 0.70221]
[-2.18677, 2.21791, 0.57243]
[-3.16291, 2.32609, 1.03614]
[-2.98868, 2.4629, 0.40026]
[-2.37152, 2.55224, 0.64066]
[-3.01113, 2.52921, 1.17777]
ratio: 20
[-4.48405, 3.00463, 1.45081]
[-2.78745, 2.25786, 0.62019]
[-3.80631, 2.71905, 1.10831]
[-2.33652, 2.0821, 0.60125]
[-2.78862, 2.57325, 0.5167]
[-1.23831, 2.17586, 0.85325]
[-3.00064, 2.318, 1.66552]
[-2.54605, 2.8181, 0.56879]
[-2.62482, 2.64067, 0.43063]
[-2.9084, 2.34565, 1.14212]
ratio: 50
[-3.61535, 2.8079, 1.90597]
[-1.765, 2.1605, 1.30307]
[-2.90974, 2.74135, 0.71409]
[-1.00802, 2.01203, 0.95008]
[-2.69431, 2.86702, 1.32219]
[-1.26446, 1.80058, 0.77125]
[-1.10304, 2.69649, 2.32091]
[-1.15043, 2.74996, 1.84605]
[-2.011, 2.66888, 0.8709]
[-2.50288, 2.23787, 1.62665]
ratio: 100
[-2.38678, 3.55555, 2.05781]
[-1.07891, 2.12761, 2.53363]
[-1.02447, 2.88643, 2.16754]
[-0.47467, 1.97916, 1.08508]
[-3.95229, 3.45432, 2.45714]
[-1.86771, 1.53591, 0.67398]
[-0.29807, 2.81919, 2.541]
[-1.16426, 2.69413, 1.97518]
[-1.24353, 2.35371, 0.90051]
[-2.73079, 2.96275, 1.71517]
ratio: 200
[-2.57214, 3.65385, 3.48845]
[-0.68745, 2.47665, 2.19183]
[2.5264, 3.1807, 5.11516]
[-0.41136, 2.27868, 0.87186]
[-2.83457, 2.97796, 3.13823]
[-2.42656, 2.67244, 1.94582]
[0.97235, 2.61841, 6.52571]
[-0.15809, 3.33408, 3.0626]
[0.37536, 4.97194, 3.57977]
[-1.90523, 2.93487, 1.37925]
ある程度は防ぐことが出来たが、切片を求めたいのでなければ、特に重みをつけるメリットはなさそう。
class-weight を使った補正の有効性はなさそう
- negative down sampling としては、重みをつけることメリットはなさそう。
balanced
を使っていいかbalanced
を logistic regression で使ってしまうと、切片重みの情報が完全に失われてしまう。
[https://swarbrickjones.wordpress.com/2017/03/28/cross-entropy-and-training-test-class-imbalance/:embed:cite]
の最後に、
Also note that we heavily relied on the assumption that the positives/negatives in X were equidistributed to the positives/negatives in X' respectively. If that it is not true, then this analysis could be of limited use!
とある通り、本来は balanced
を使うべきではないだろう。
ところが、down sampling できる状況では、切片重みは事前に求められるはずである。むしろ、ダウンサンプリングして切片重み求めることのほうが効率が悪いかもしれない。切片重みは全体の事例から求める( [tex:{ \ln \hat{pi} }] )として、balanced
オプションで、切片は0としてfit_intercept=False
に設定し計算してみる。
def check_coefs_for_each_ratio_with_balanced(df_list): for ratio in [1, 2, 5, 10, 20, 50, 100, 200]: display(Markdown("### ratio: %s" % str(ratio))) for i in range(0, 10): df = df_list[i] intercept = np.log( float(len(df[df['y'] == 1])) / len(df) ) df_d = df[df['r'] < 1.0/ratio] logreg = LogisticRegression(penalty='l1', C=1, fit_intercept=False, class_weight='balanced') column_list = ['x'] + ['z'+str(i) for i in range(FAKE)] logreg.fit(df_d[column_list], pd.Series(df_d['y'].values.flatten())) w_x = logreg.coef_.tolist()[0][0] w_z = sum([w*w for w in logreg.coef_.tolist()[0][1:]]) # 二乗して合計(ここはフェイクなので小さい方がいい) print([round(x,5) for x in [intercept, w_x, w_z] ])
check_coefs_for_each_ratio_with_balanced(df_list)
ratio: 1
[-4.19971, 2.80006, 0.58413]
[-4.14144, 2.37942, 1.11962]
[-4.14775, 2.60704, 0.80826]
[-4.05129, 2.12227, 0.7282]
[-4.43966, 2.70717, 1.61314]
[-4.2475, 2.38503, 0.50004]
[-4.23361, 2.15592, 1.2708]
[-4.2687, 2.61576, 0.52414]
[-4.16048, 2.54132, 0.72488]
[-4.20639, 2.43589, 1.23491]
ratio: 2
[-4.19971, 2.73953, 0.51565]
[-4.14144, 2.32173, 1.09193]
[-4.14775, 2.59147, 0.94536]
[-4.05129, 2.08336, 0.76821]
[-4.43966, 2.6852, 1.72355]
[-4.2475, 2.35151, 0.54959]
[-4.23361, 2.19345, 1.20504]
[-4.2687, 2.57608, 0.56743]
[-4.16048, 2.48444, 0.72722]
[-4.20639, 2.41925, 1.15235]
ratio: 5
[-4.19971, 2.73559, 0.41098]
[-4.14144, 2.40825, 0.89882]
[-4.14775, 2.51915, 1.05386]
[-4.05129, 2.29807, 0.76]
[-4.43966, 2.63634, 1.39499]
[-4.2475, 2.27439, 0.62931]
[-4.23361, 2.21276, 1.2259]
[-4.2687, 2.58391, 0.63055]
[-4.16048, 2.61542, 0.69615]
[-4.20639, 2.27213, 1.08081]
ratio: 10
[-4.19971, 2.72486, 0.48284]
[-4.14144, 2.38394, 0.65321]
[-4.14775, 2.47043, 1.16771]
[-4.05129, 2.2557, 0.76595]
[-4.43966, 2.64495, 1.43846]
[-4.2475, 2.17495, 0.47283]
[-4.23361, 2.31494, 1.26874]
[-4.2687, 2.49734, 0.62739]
[-4.16048, 2.55463, 0.55818]
[-4.20639, 2.50076, 1.1878]
ratio: 20
[-4.19971, 2.79394, 0.49988]
[-4.14144, 2.20952, 0.49212]
[-4.14775, 2.61027, 0.60249]
[-4.05129, 2.04884, 0.57751]
[-4.43966, 2.53871, 0.77375]
[-4.2475, 2.09275, 0.34751]
[-4.23361, 2.24294, 1.47811]
[-4.2687, 2.77854, 0.68346]
[-4.16048, 2.61743, 0.496]
[-4.20639, 2.29761, 1.10332]
ratio: 50
[-4.19971, 2.57246, 0.9157]
[-4.14144, 2.07143, 0.86454]
[-4.14775, 2.63409, 0.61976]
[-4.05129, 1.91869, 0.65903]
[-4.43966, 2.65093, 0.9455]
[-4.2475, 1.69904, 0.30351]
[-4.23361, 2.48733, 1.46496]
[-4.2687, 2.55024, 0.71679]
[-4.16048, 2.56966, 0.81352]
[-4.20639, 2.03153, 0.87826]
ratio: 100
[-4.19971, 3.18644, 0.98824]
[-4.14144, 1.96485, 1.55056]
[-4.14775, 2.66571, 1.46552]
[-4.05129, 1.83133, 0.55313]
[-4.43966, 2.9996, 1.15968]
[-4.2475, 1.3735, 0.21852]
[-4.23361, 2.57707, 1.42351]
[-4.2687, 2.45839, 0.84825]
[-4.16048, 2.25478, 0.61956]
[-4.20639, 2.65479, 0.88477]
ratio: 200
[-4.19971, 3.15883, 1.66064]
[-4.14144, 2.32423, 1.09636]
[-4.14775, 2.85993, 2.00946]
[-4.05129, 2.1121, 0.3894]
[-4.43966, 2.5296, 1.53771]
[-4.2475, 2.26719, 0.72349]
[-4.23361, 2.36902, 2.97373]
[-4.2687, 3.02366, 1.29845]
[-4.16048, 3.97411, 1.76377]
[-4.20639, 2.64598, 0.87854]
よさそうである。まず切片はダウンサンプリングに関係なく全体事例を使っているのでかなり精度よく推定されているといえる。x
の重み(2.302585が理想)もこれまでの方法の中で最も良い値が得られた様に思える。これは切片が固定されているため、これが探索空間の絞り込みになり非常に良い推定を与えると考えられる。
balanaced
の学習パラメータの探索
def check_learn_param_balanced(): for ratio in [1, 2, 5, 10, 20, 50, 100, 200]: display(Markdown("### ratio: %s" % str(ratio))) df_d = df[df['r'] < 1.0/ratio] column_list = ['x'] + ['z'+str(i) for i in range(FAKE)] grid_search_result = grid_search_cv(df_d[column_list], pd.Series(df_d['y'].values.flatten()), class_weight='balanced', fit_intercept=False) df_result = pd.DataFrame(grid_search_result.cv_results_) display(df_result[['params', 'mean_test_neg_log_loss', 'std_test_neg_log_loss']]) check_learn_param_balanced()
ratio: 1
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.579756 | 0.016201 |
1 | {u'C': 0.1} | -0.557804 | 0.020383 |
2 | {u'C': 1} | -0.557315 | 0.020976 |
3 | {u'C': 10} | -0.557317 | 0.021046 |
4 | {u'C': 100} | -0.557318 | 0.021053 |
ratio: 2
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.614889 | 0.015515 |
1 | {u'C': 0.1} | -0.562747 | 0.024808 |
2 | {u'C': 1} | -0.561251 | 0.026736 |
3 | {u'C': 10} | -0.561220 | 0.026921 |
4 | {u'C': 100} | -0.561219 | 0.026940 |
ratio: 5
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.704480 | 0.007989 |
1 | {u'C': 0.1} | -0.569775 | 0.027439 |
2 | {u'C': 1} | -0.565177 | 0.031148 |
3 | {u'C': 10} | -0.565230 | 0.031662 |
4 | {u'C': 100} | -0.565249 | 0.031717 |
ratio: 10
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.692749 | 0.000850 |
1 | {u'C': 0.1} | -0.582792 | 0.040261 |
2 | {u'C': 1} | -0.570911 | 0.051377 |
3 | {u'C': 10} | -0.571530 | 0.053412 |
4 | {u'C': 100} | -0.571638 | 0.053624 |
ratio: 20
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.693147 | 2.135512e-16 |
1 | {u'C': 0.1} | -0.591359 | 5.392871e-02 |
2 | {u'C': 1} | -0.573781 | 8.087916e-02 |
3 | {u'C': 10} | -0.575191 | 8.687242e-02 |
4 | {u'C': 100} | -0.575433 | 8.755660e-02 |
ratio: 50
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.693147 | 1.110223e-16 |
1 | {u'C': 0.1} | -0.610073 | 4.553321e-02 |
2 | {u'C': 1} | -0.589434 | 7.575629e-02 |
3 | {u'C': 10} | -0.598316 | 8.391099e-02 |
4 | {u'C': 100} | -0.599609 | 8.501846e-02 |
ratio: 100
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.693147 | 2.220446e-16 |
1 | {u'C': 0.1} | -0.588712 | 4.332537e-02 |
2 | {u'C': 1} | -0.574640 | 1.232090e-01 |
3 | {u'C': 10} | -0.598204 | 1.625511e-01 |
4 | {u'C': 100} | -0.601811 | 1.678043e-01 |
ratio: 200
params | mean_test_neg_log_loss | std_test_neg_log_loss | |
---|---|---|---|
0 | {u'C': 0.01} | -0.693147 | 1.110223e-16 |
1 | {u'C': 0.1} | -0.580985 | 3.726425e-02 |
2 | {u'C': 1} | -0.583421 | 1.358876e-01 |
3 | {u'C': 10} | -0.626849 | 2.029742e-01 |
4 | {u'C': 100} | -0.643918 | 2.217282e-01 |
特に極値となる学習パラメータが変動することもなさそうで、よさそう。
Logistic Regression における Negative Down Sampling の設定まとめ
balanced
を使う方法- 切片重みを down sampling する前に求める
class_weight='balanced', fit_intercept=False
で学習を行うRe-Calibration
は必要ない
- どうしても
balanced
を使いたくない場合intercept_scale=10000
にするRe-Calibration
が必要
お約束
小標本における正規分布の信頼区間を Excel で求める
背景
コイン投げのようなベルヌーイ試行の確率の信頼区間は
Wilson score interval を使う。 - 中野智文のブログ
でいいとして、母集団が正規分布で小標本の場合は、 t 分布で求めるのが、一般的である。
母集団が正規分布に従うときで標本の大きさ(サンプルサイズ)が大きいときは、Wald法による信頼区間を使うことができる。 その標本の大きさが小さいときWald法は使えない。具体的には30未満のときはWald法は使えない。このような Wald法が使えない標本の大きさを 小標本という。 (ベルヌーイ分布の場合には、標本平均にもよるが、もっと大きな数が必要) 母集団が正規分布に従い、小標本だった場合には、t分布を用いる必要がある。
Excel で t 分布 の信頼区間
Excelでは、CONFIDENCE.T 関数 - Office サポート という関数により、t 分布で信頼区間を求めることができる。
一見そのまま使えば良いようだが、標準偏差とは、不偏標準偏差のことだろうか、それとも標本標準偏差のことだろうか。
答え
不偏標準偏差 を使う。
確認
確認は、確率統計のはなし の「貯金額の標本が10万円と30万円の二人だった場合、全員の平均値の信頼区間」の計算結果 p.47
μ = 20±63.14
にて確認する。(なおこの本では、標本標準偏差を使っているが、気にしない)
今、B1 と、B2 にそれぞれ10、30 が入っている。
=CONFIDENCE.T(0.1,STDEV(B1:B2),COUNT(B1:B2))
を計算すると、63.13751515が得られるはず。
ここで、STDEV
は不偏標準偏差を、COUNT
は標本の大きさを示す。
Google SpreadSheet では?
CONFIDENCE.T
がない。CONFIDENCE
に統合されているのかと思ったらそうでもないみたい。
信頼区間は次の式で求まる。
ここで、 は、t分布表の t の値、 は標本標準偏差を表す。
Google SpreadSheet では、t 分布表は、
TINV - ドキュメント エディタ ヘルプ から得る。例えば、先の例だと、10%で、標本の大きさは2なので、自由度はマイナス1して、1。
よって、TINV(0.1, 1)
で、6.313751515と計算されるはず。
Google SpreadSheet では、標本標準偏差は、
STDEVP - ドキュメント エディタ ヘルプ から得る。例えば、先の例だと、STDEVP(B1:B2)
で10.0と計算されるはず。
結局
=TINV(0.1,COUNT(B1:B2)-1)*STDEVP(B1:B2)
となり、63.13751515 と表示されれば、確認完了。
念のため python での求め方
from scipy import stats from scipy.stats import t import numpy as np
a = [10.0, 30.0]
t0, t1 = stats.t.interval(alpha=0.90, df=len(a)-1, loc=np.mean(a), scale=stats.sem(a)) print([t0, t1])
[-43.137515148009399, 83.137515148009328]
信頼区間を直接求めるので、こうなる。
補足
ところで、Microsoft のサイトでは、標本数とあるが、標本の大きさの間違いだろう
なお、確率統計のはなし では標本の大きさとして 標本の数 というのが出てくるが、「の」が間に入っているからセーフ?
scikit learn で cross validation で confusion matrix を取得する。
背景
confusion matrix を取得する場合は、一部の例だけでなく cross validation で全ての事例に対して取得したい。
対応
sklearn.model_selection.cross_val_predict を使う。
例
今回は iris のデータを使う。当たり前だが、confusion matrix は普通分類問題に使うよね?
from sklearn import datasets iris = datasets.load_iris() X = iris.data y = iris.target
にも関わらず、例だから、logistic regression などを使ってみる。
from sklearn.linear_model import LogisticRegression logistic = LogisticRegression()
ここで、cross validation 登場。 cv
パラメータは、KFold(n_splits=10, shuffle=True)
に設定している*1。推定結果と教師データから、confusion matrix を生成する。
from sklearn.model_selection import (cross_val_predict, KFold) from sklearn.metrics import confusion_matrix y_pred = cross_val_predict(logistic,X,y,cv=KFold(n_splits=10, shuffle=True)) conf_mat = confusion_matrix(y,y_pred)
コードを書かなくていいし、そのコードのバグも心配しなくていい*2。検証するコードがバグっていて、何を色々努力してもダメ、というのは最悪なパターンだと思うけど、 機械学習ではよくある話のような気がする。
そして、よくある、 confusion matrix を plot するコード。
# scikit-learn.org/stable/auto_examples/model_selection/plot_confusion_matrix.html import matplotlib.pyplot as plt %matplotlib inline import itertools import numpy as np def plot_confusion_matrix(cm, classes, normalize=False, title='Confusion matrix', cmap=plt.cm.Blues): """ This function prints and plots the confusion matrix. Normalization can be applied by setting `normalize=True`. """ if normalize: cm = cm.astype('float') / cm.sum(axis=1)[:, np.newaxis] print("Normalized confusion matrix") else: print('Confusion matrix, without normalization') print(cm) plt.imshow(cm, interpolation='nearest', cmap=cmap) plt.title(title) plt.colorbar() tick_marks = np.arange(len(classes)) plt.xticks(tick_marks, classes, rotation=45) plt.yticks(tick_marks, classes) fmt = '.2f' if normalize else 'd' thresh = cm.max() / 2. for i, j in itertools.product(range(cm.shape[0]), range(cm.shape[1])): plt.text(j, i, format(cm[i, j], fmt), horizontalalignment="center", color="white" if cm[i, j] > thresh else "black") plt.tight_layout() plt.ylabel('True label') plt.xlabel('Predicted label')
そして実行+結果
plot_confusion_matrix(conf_mat, list(iris.target_names))
Confusion matrix, without normalization
[[50 0 0]
[ 0 45 5]
[ 0 2 48]]
まとめ
confusion matrix を取得したいときは、sklearn.model_selection.cross_val_predict を使う。 このcolab
広告
*1:詳しくは scikit learn の Kfold, StratifiedKFold, ShuffleSplit の違い - 中野智文のブログ
*2:と書いたけどこれは嘘。心配(注意)はした方がいい。ただ有名なコードなら多くの人が心配してくれるし、直してもくれるから安心
scikit learn で DataFrameから色々なスパースな型に変換して学習
背景
scikit learn で学習しようとすると、メモリーを使い尽くす。
色々なスパースな型に変換して学習
準備
まず、データは次のものを利用。
measurements = [ {'city': 'Dubai', 'temperature': 31.0, 'country': 'U.A.E.'}, {'city': 'London', 'country': 'U.K.', 'temperature': 27.0}, {'city': 'San Fransisco', 'country': 'U.S.', 'temperature': 24.0}, ]
ただし変換をして、DataFrame用に。
x_column_dict = {} for d in measurements: for c, v in d.items(): if c in x_column_dict: x_column_dict[c].append(v) else: x_column_dict[c] = [v] x_column_dict
{'city': ['Dubai', 'London', 'San Fransisco'],
'country': ['U.A.E.', 'U.K.', 'U.S.'],
'temperature': [31.0, 27.0, 24.0]}
import pandas as pd df = pd.DataFrame(x_column_dict) df
city | country | temperature | |
---|---|---|---|
0 | Dubai | U.A.E. | 31.0 |
1 | London | U.K. | 27.0 |
2 | San Fransisco | U.S. | 24.0 |
教師データを準備
import numpy as np y = np.array([1,0,0])
get_dummies で sparse=True
を使う。
X = pd.get_dummies(df, sparse=True)
X.info()
<class 'pandas.core.frame.DataFrame'>
RangeIndex: 3 entries, 0 to 2
Data columns (total 7 columns):
temperature 3 non-null float64
city_Dubai 3 non-null uint8
city_London 3 non-null uint8
city_San Fransisco 3 non-null uint8
country_U.A.E. 3 non-null uint8
country_U.K. 3 non-null uint8
country_U.S. 3 non-null uint8
dtypes: float64(1), uint8(6)
memory usage: 114.0 bytes
↑ SparseDataFrame
型ではないのですが、詐欺ですか?
from sklearn.linear_model import LogisticRegression logistic = LogisticRegression() logistic.fit(X, y)
LogisticRegression(C=1.0, class_weight=None, dual=False, fit_intercept=True,
intercept_scaling=1, max_iter=100, multi_class='ovr', n_jobs=1,
penalty='l2', random_state=None, solver='liblinear', tol=0.0001,
verbose=0, warm_start=False)
一応学習はできる(DataFrame
型なので)
to_sparse()
を使う。
X = pd.get_dummies(df).to_sparse() X.info()
<class 'pandas.core.sparse.frame.SparseDataFrame'>
RangeIndex: 3 entries, 0 to 2
Data columns (total 7 columns):
temperature 3 non-null float64
city_Dubai 3 non-null uint8
city_London 3 non-null uint8
city_San Fransisco 3 non-null uint8
country_U.A.E. 3 non-null uint8
country_U.K. 3 non-null uint8
country_U.S. 3 non-null uint8
dtypes: float64(1), uint8(6)
memory usage: 114.0 bytes
確かに SparseDataFrame
型に変換されている。ただしメモリの使用量は同じ。
logistic.fit(X, y)
LogisticRegression(C=1.0, class_weight=None, dual=False, fit_intercept=True,
intercept_scaling=1, max_iter=100, multi_class='ovr', n_jobs=1,
penalty='l2', random_state=None, solver='liblinear', tol=0.0001,
verbose=0, warm_start=False)
学習できたことを確認。
さらに、fill_value=0
を使う。
X = pd.get_dummies(df).to_sparse(fill_value=0)
X.info()
<class 'pandas.core.sparse.frame.SparseDataFrame'>
RangeIndex: 3 entries, 0 to 2
Data columns (total 7 columns):
temperature 3 non-null float64
city_Dubai 3 non-null uint8
city_London 3 non-null uint8
city_San Fransisco 3 non-null uint8
country_U.A.E. 3 non-null uint8
country_U.K. 3 non-null uint8
country_U.S. 3 non-null uint8
dtypes: float64(1), uint8(6)
memory usage: 102.0 bytes
↑メモリの使用量が減っている。ただしこれは学習前の話。
logistic.fit(X, y)
LogisticRegression(C=1.0, class_weight=None, dual=False, fit_intercept=True,
intercept_scaling=1, max_iter=100, multi_class='ovr', n_jobs=1,
penalty='l2', random_state=None, solver='liblinear', tol=0.0001,
verbose=0, warm_start=False)
↑実際には上記の処理を行っても、かなりメモリを使ってしまう。
New!! さらに coo に変換する。
[https://pandas.pydata.org/pandas-docs/stable/generated/pandas.SparseSeries.to_coo.html]
さらに、scipy の sparse matrix に変換できるらしい。(@hagino3000 さんに教えてもらった)
X = pd.get_dummies(df).to_sparse(fill_value=0) Xcoo = X.to_coo() print Xcoo
(0, 0) 31.0
(1, 0) 27.0
(2, 0) 24.0
(0, 1) 1.0
(1, 2) 1.0
(2, 3) 1.0
(0, 4) 1.0
(1, 5) 1.0
(2, 6) 1.0
logistic.fit(Xcoo, y)
LogisticRegression(C=1.0, class_weight=None, dual=False, fit_intercept=True,
intercept_scaling=1, max_iter=100, multi_class='ovr', n_jobs=1,
penalty='l2', random_state=None, solver='liblinear', tol=0.0001,
verbose=0, warm_start=False)
↑学習時のメモリの使用量がぐっと減る。
つづき
logistic.coef_
array([[-0.01603272, 0.43215653, -0.25978198, -0.26661549, 0.43215653,
-0.25978198, -0.26661549]])
w = dict(zip(X.columns.tolist(), logistic.coef_[0])) w
{'city_Dubai': 0.43215653143933291,
'city_London': -0.25978198285599019,
'city_San Fransisco': -0.26661548568954441,
'country_U.A.E.': 0.43215653143933291,
'country_U.K.': -0.25978198285599019,
'country_U.S.': -0.26661548568954441,
'temperature': -0.016032719041480594}
b = logistic.intercept_[0]
b
-0.094240937106201669
まとめ
get_dummies
のsparse=True
オプションは意味がないto_sparse(fill_value=0)
で0を除去。to_coo()
で学習時のメモリを最小化。(@hagino3000 さん感謝)