001package org.opengion.penguin.math.statistics; 002 003import java.util.Arrays; 004import java.util.Collections; 005import java.util.List; 006 007/** 008 * 多項ロジスティック回帰の実装です。 009 * 確率的勾配降下法(SGD)を利用します。 010 * 011 * ロジスティック回帰はn次元の情報からどのグループに所属するかの予測値を得るための手法の一つです。 012 * 013 * 実装は 014 * http://nbviewer.jupyter.org/gist/mitmul/9283713 015 * https://yusugomori.com/projects/deep-learning/ 016 * を参考にしています。 017 */ 018public class HybsLogisticRegression { 019 private int n_N; // データ個数 020 private int n_in; // データ次元 021 private int n_out; // ラベル種別数 022 private Integer[] random_index; //確率的勾配降下用index 023 024 // 写像変数ベクトル f(x) = Wx + b 025 private double[][] vW; 026 private double[] vb; 027 028 029 030 /** 031 * コンストラクタ。 032 * 033 * 学習もしてしまう。 034 * 035 * xはデータセット各行がn次元の説明変数となっている。 036 * trainはそれに対する{0,1,0},{1,0,0}のようなラベルを示すベクトルとなる。 037 * 学習率は通常、0.1程度を設定する。 038 * このロジックではループ毎に0.95をかけて徐々に学習率が下がるようにしている。 039 * 全データを利用すると時間がかかる場合があるので、確率的勾配降下法を利用しているが、 040 * 選択個数はデータに対する割合を与える。 041 * データ個数が少ない場合は1をセットすればよい。 042 * 043 * @param data データセット配列 044 * @param label データに対応したラベルを示す配列 045 * @param learning_rate 学習係数(0から1の間の数値) 046 * @param loop 学習のループ回数(ミニバッチを作る回数) 047 * @param minibatch_rate 全体に対するミニバッチの割合(0から1の間の数値) 048 * 049 */ 050 public HybsLogisticRegression(double data[][], int label[][], double learning_rate ,int loop, double minibatch_rate ) { 051 List<Integer> indexList; //シャッフル用 052 053 this.n_N = data.length; 054 this.n_in = data[0].length; 055 this.n_out = label[0].length; // ラベル種別 056 057 vW = new double[n_out][n_in]; 058 vb = new double[n_out]; 059 060 // 確率勾配に利用するための配列インデックス配列 061 random_index = new Integer[n_N]; //プリミティブ型だとasListできないため 062 for( int i=0; i<n_N; i++){ 063 random_index[i] = i; 064 } 065 indexList = Arrays.asList( random_index ); 066 067 for(int epoch=0; epoch<loop; epoch++) { 068 Collections.shuffle( indexList ); 069 random_index = (Integer[])indexList.toArray(new Integer[indexList.size()]); 070 071 //random_indexの先頭からn_N*minibatch_rate個のものを対象に学習をかける(ミニバッチ) 072 for(int i=0; i< n_N * minibatch_rate; i++) { 073 int idx = random_index[i]; 074 train(data[idx], label[idx], learning_rate); 075 } 076 learning_rate *= 0.95; //徐々に学習率を下げて振動を抑える。 077 } 078 } 079 080 /** 081 * コンストラクタ。 082 * 083 * Wとbのセットのみを行い、過去の結果から予測値計算だけを行う場合。 084 * 085 * @param W_in 係数 086 * @param b_in バイアス 087 * 088 */ 089 public HybsLogisticRegression(double W_in[][], double[] b_in) { 090 this.n_in = W_in[0].length; 091 this.n_out = b_in.length; // ラベル種別 092 093 // ディープコピーはしない 094 vW = W_in; 095 vb = b_in; 096 } 097 098 099 /** 100 * データを与えて学習をさせます。 101 * パラメータの1行を与えています。 102 * 103 * 0/1のロジスティック回帰の場合は 104 * ラベルc(0or1)が各xに対して与えられている時 105 * s(x)=σ(Wx+b)=1/(1+ exp(-Wx-b))として、 106 * 確率の対数和L(W,b)の符号反転させたものの偏導関数 107 * ∂L/∂w=-肺(c-s(x)) 108 * ∂L/∂b=-=(c-s(x)) 109 * が最小になるようなW,bの値をパラメータを変えながら求める。 110 * というのが実装になる。(=0を求められないため) 111 * 多次元の場合はシグモイド関数σ(x)の代わりにソフトマックス関数π(x)を利用して 112 * 拡張したものとなる。(以下はソフトマックス関数利用) 113 * 114 * @param in_x 1行分のデータ 115 * @param in_y xに対するラベル 116 * @param lr 学習率 117 * @return 差分配列 118 */ 119 private double[] train(double[] in_x, int[] in_y, double lr) { 120 double[] p_y_given_x = new double[n_out]; 121 double[] dy = new double[n_out]; 122 123 for(int i=0; i<n_out; i++) { 124 p_y_given_x[i] = 0; 125 for(int j=0; j<n_in; j++) { 126 p_y_given_x[i] += vW[i][j] * in_x[j]; 127 } 128 p_y_given_x[i] += vb[i]; 129 } 130 softmax(p_y_given_x); 131 132 // 勾配の平均で更新? 133 for(int i=0; i<n_out; i++) { 134 dy[i] = in_y[i] - p_y_given_x[i]; 135 136 for(int j=0; j<n_in; j++) { 137 vW[i][j] += lr * dy[i] * in_x[j] / n_N; 138 } 139 140 vb[i] += lr * dy[i] / n_N; 141 } 142 143 return dy; 144 } 145 146 /** 147 * ソフトマックス関数 148 * π(xi) = exp(xi)/Σexp(x) 149 * @param in_x 変数X 150 */ 151 private void softmax(double[] in_x) { 152 // double max = 0.0; 153 double sum = 0.0; 154 155 // for(int i=0; i<n_out; i++){ 156 // if(max < x[i]){ 157 // max = x[i]; 158 // } 159 // } 160 161 for(int i=0; i<n_out; i++) { 162 //x[i] = Math.exp(x[i] - max); // maxとの差分を取ると利点があるのか分からなかった 163 in_x[i] = Math.exp(in_x[i]); 164 sum += in_x[i]; 165 } 166 167 for(int i=0; i<n_out; i++){ 168 in_x[i] /= sum; 169 } 170 } 171 172 /** 173 * 写像式 Wx+b のW、係数ベクトル。 174 * @return 係数ベクトル 175 */ 176 public double[][] getW(){ 177 return vW; 178 } 179 180 /** 181 * 写像式 Wx + bのb、バイアス。 182 * @return バイアスベクトル 183 */ 184 public double[] getB(){ 185 return vb; 186 } 187 188 /** 189 * 出来た予測式に対して、データを入力してyを出力する。 190 * (yは各ラベルに対する確率分布となる) 191 * @param in_x 予測したいデータ 192 * @return 予測結果 193 */ 194 public double[] predict(double[] in_x) { 195 double[] out_y = new double[n_out]; 196 197 for(int i=0; i<n_out; i++) { 198 out_y[i] = 0.; 199 for(int j=0; j<n_in; j++) { 200 out_y[i] += vW[i][j] * in_x[j]; 201 } 202 out_y[i] += vb[i]; 203 } 204 205 softmax(out_y); 206 207 return out_y; 208 } 209 210 /*** ここまでが本体 ***/ 211 /*** ここからテスト用mainメソッド ***/ 212 /** 213 * @param args 214 *****************************************/ 215 public static void main(String[] args) { 216 // 3つの分類で分ける 217 double[][] train_X = { 218 {-2.0, 2.0} 219 ,{-2.1, 1.9} 220 ,{-1.8, 2.1} 221 ,{0.0, 0.0} 222 ,{0.2, -0.2} 223 ,{-0.1, 0.1} 224 ,{2.0, -2.0} 225 ,{2.2, -2.1} 226 ,{1.9, -2.0} 227 }; 228 229 int[][] train_Y = { 230 {1, 0, 0} 231 ,{1, 0, 0} 232 ,{1, 0, 0} 233 ,{0, 1, 0} 234 ,{0, 1, 0} 235 ,{0, 1, 0} 236 ,{0, 0, 1} 237 ,{0, 0, 1} 238 ,{0, 0, 1} 239 }; 240 241 // test data 242 double[][] test_X = { 243 {-2.5, 2.0} 244 ,{0.1, -0.1} 245 ,{1.5,-2.5} 246 }; 247 248 double[][] test_Y = new double[test_X.length][train_Y[0].length]; 249 250 HybsLogisticRegression hlr = new HybsLogisticRegression( train_X, train_Y, 0.1, 500, 1 ); 251 252 // テスト 253 // このデータでは2番目の条件には入りにくい? 254 for(int i=0; i<test_X.length; i++) { 255 test_Y[i] = hlr.predict(test_X[i]); 256 System.out.print( Arrays.toString(test_Y[i]) ); 257 } 258 } 259} 260