xgboost: Higgs Boson Machine Learning Challenge
13474 ワード
コードの元の出典:https://github.com/dmlc/xgboost/tree/master/demo/kaggle-higgs
一、問題の紹介
コンテスト公式サイト:https://www.kaggle.com/c/higgs-boson/
ヒッグスボクロトン(英語:Higgs boson)は標準モデルの基本粒子で、物理学者ピーター・ヒッグスによって命名された.
2012年7月4日、欧州原子力研究機構(CERN)は、LHCのコンパクトな微小子コイル(CMS)が125.3±0.6 GeVの新ボクロトン(バックグラウンドの期待値4.9の標準差を超える)を検出し、超環状機器(ATLAS)が126.5 GeVの新ボクロトン(5の標準差)を測定し、この2つの粒子はヒッグスボクロに極めて似ていると発表した.
2013年3月14日、欧州の核研究機関はニュース原稿を発表し、以前に検出された新しい粒子は一時的にヒグスボクロトンであることが確認され、ゼロスピンと双宇を持っていると発表した.これはヒグスボクロトンが持つべき2つの基本的な性質だが、一部の実験の結果は理論予測に合致せず、より多くのデータは処理と分析を待っている.2013年10月8日、「二次原子粒子質量の生成メカニズム理論は、人類のこの方面に対する理解を促進し、最近、ヨーロッパの核研究組織傘下の大型強子衝突機の超環面機器とコンパクトμサブコイル検出器で発見された基本粒子が確認された」と、フランソワ・エングラー、ピーター・ヒッグスが2013年のノーベル物理学賞を受賞した.1つの粒子の重要な特徴は、他の粒子の後にどれだけ遅延するかである.CERNはATLASを用いて物理実験を行い新しい粒子を探す.実験は最近,二つのtau粒子でHiggs boson遅延が現れることを見出したが,この遅延は背景雑音に埋もれた小さな信号にすぎなかった.
このコンテストの目的は,機械学習法を用いて,ATLAS実験による粒子発見の顕著性を高めることである.競合には粒子物理の背景知識は必要ありません(実際の問題を解決する際に背景知識は大きく役立ちます).コンテストデータは、ATLASが検出したイベントの特徴に基づいて合成されたデータであり、コンテストタスクは、イベントを「tau tau decay of a Higgs boson」または「background」に分類することである
二、データセットの概要
トレーニングデータセット:http://download.csdn.net/download/u011630575/10268125
training.csv:トレーニングセットには250000個のイベントが含まれており、各イベントにはID、30個の特徴、重み、ラベルがあります.
test.csv:テストデータには550000イベントが含まれており、各イベントにはIDと30の特徴が含まれています.すべての変数はPRI_を除いてfloating pointタイプです.jet_numはintegerがPRI(PRImitives)を用いたプレフィックス特徴を検出器として測定するbunch collision「オリジナル」に関するデータである.DER(DERived)をATLASとする物理学者が選択した元の特徴から計算したデータ.欠落データは−999.0と表記され,すべての特徴の正常値とは異なる.
特徴、重み、およびラベルの特定の意味は、CERNの技術文書を参照することができる.(コンテスト公式サイトにリンクあり)
三、コード説明
1、データの読み取り
2、データの前処理
3、訓練環境の準備
トレーニングパラメータの設定:1.objective[デフォルトreg:linear]:最小化する必要がある損失関数を定義します:binary:logistic二分類の論理回帰、予測を返す確率(カテゴリではありません).Higgs Bosonコンテストは2種類の分類タスクであり,2分類の論理回帰を採用している.2. scale_pos_Weight[デフォルト1]:各種類の別サンプルが非常にアンバランスである場合、このパラメータを正の値に設定すると、アルゴリズムをより速く収束させることができる.Higgs Bosonコンテストでは,各(正/負)サンプルの重みが与えられ,すべての正/負サンプルの重みが加算され,トレーニングセットの正負サンプルの割合が得られる.3.eta[デフォルト0.3]は学習率です.オーバーフィットを防止するために、更新中に使用される収縮ステップ長.計算を上げるたびに、アルゴリズムは新しいフィーチャーの重みを直接取得します.etaは,特徴の重みを縮小することによって,リフト計算過程をより保守的にする.値の範囲は[0,1]4である.max_depth[デフォルト6]は、ツリーの最大深さを定義します.この値は、オーバーフィットを回避するためにも使用されます.max_depthが大きいほどモデルが複雑になり,より具体的で局所的なサンプルを学ぶ.標準:3-10 5.eval_metric[デフォルト値はobjectiveパラメータの値に依存する]有効なデータのメトリックメソッド.回帰問題の場合、デフォルト値はrmse、分類問題の場合、デフォルト値はerrorです.典型的な値は、rmse平均二乗誤差mae平均絶対誤差logloss負対数尤度関数値error二分類誤り率(閾値0.5)merror多分類誤り率mlogloss多分類logloss損失関数auc曲線下面積(Area Under Curve):異なるしきい値でのモデルの正解率である.6.nthread[デフォルト値は最大可能なスレッド数]このパラメータはマルチスレッド制御に使用され、システムのコア数を入力する必要があります.
CPUのすべてのコアを使用したい場合は、このパラメータを入力しないでください.アルゴリズムは自動的に検出します.
トレーニングデータはDMatixにインポートされ、後続のトレーニングがより速くなります.
4、xgboostモデルトレーニング、保存
train関数のパラメータ:xgboost.train(params,dtrain,num_boost_round=10,evals=(),obj=None,feval=None,maximize=False,early_stopping_rounds=None, evals_result=None,verbose_eval=True,learning_rates=None,xgb_model=None)
params:これは辞書で、訓練中のパラメータキーワードと対応する値が含まれています.形式はparams={‘booster’:’gbtree’,’eta’:0.1}です.
dtrainトレーニングのデータ
num_boost_roundこれはリフト反復の個数を指します
evals:トレーニング中にリストの要素を評価するためのリストです.形式はevals=[(dtrain,’train’),(dval,’val’)]またはevals=[(dtrain,’train’)]である.
このコードでは、トレーニング中に検証セットの効果を観察できるようにする最初のケースを使用します.
5、テストプロセス
テストデータの読み取りと訓練されたモデル
XGBoost試験環境準備
テストデータインポートDMatrixモデルインポート
テスト
テスト結果を整理し、結果を書き込んでファイルを提出する
四、コード整理
モデル構築コード
テストコード
一、問題の紹介
コンテスト公式サイト:https://www.kaggle.com/c/higgs-boson/
ヒッグスボクロトン(英語:Higgs boson)は標準モデルの基本粒子で、物理学者ピーター・ヒッグスによって命名された.
2012年7月4日、欧州原子力研究機構(CERN)は、LHCのコンパクトな微小子コイル(CMS)が125.3±0.6 GeVの新ボクロトン(バックグラウンドの期待値4.9の標準差を超える)を検出し、超環状機器(ATLAS)が126.5 GeVの新ボクロトン(5の標準差)を測定し、この2つの粒子はヒッグスボクロに極めて似ていると発表した.
2013年3月14日、欧州の核研究機関はニュース原稿を発表し、以前に検出された新しい粒子は一時的にヒグスボクロトンであることが確認され、ゼロスピンと双宇を持っていると発表した.これはヒグスボクロトンが持つべき2つの基本的な性質だが、一部の実験の結果は理論予測に合致せず、より多くのデータは処理と分析を待っている.2013年10月8日、「二次原子粒子質量の生成メカニズム理論は、人類のこの方面に対する理解を促進し、最近、ヨーロッパの核研究組織傘下の大型強子衝突機の超環面機器とコンパクトμサブコイル検出器で発見された基本粒子が確認された」と、フランソワ・エングラー、ピーター・ヒッグスが2013年のノーベル物理学賞を受賞した.1つの粒子の重要な特徴は、他の粒子の後にどれだけ遅延するかである.CERNはATLASを用いて物理実験を行い新しい粒子を探す.実験は最近,二つのtau粒子でHiggs boson遅延が現れることを見出したが,この遅延は背景雑音に埋もれた小さな信号にすぎなかった.
このコンテストの目的は,機械学習法を用いて,ATLAS実験による粒子発見の顕著性を高めることである.競合には粒子物理の背景知識は必要ありません(実際の問題を解決する際に背景知識は大きく役立ちます).コンテストデータは、ATLASが検出したイベントの特徴に基づいて合成されたデータであり、コンテストタスクは、イベントを「tau tau decay of a Higgs boson」または「background」に分類することである
二、データセットの概要
トレーニングデータセット:http://download.csdn.net/download/u011630575/10268125
training.csv:トレーニングセットには250000個のイベントが含まれており、各イベントにはID、30個の特徴、重み、ラベルがあります.
test.csv:テストデータには550000イベントが含まれており、各イベントにはIDと30の特徴が含まれています.すべての変数はPRI_を除いてfloating pointタイプです.jet_numはintegerがPRI(PRImitives)を用いたプレフィックス特徴を検出器として測定するbunch collision「オリジナル」に関するデータである.DER(DERived)をATLASとする物理学者が選択した元の特徴から計算したデータ.欠落データは−999.0と表記され,すべての特徴の正常値とは異なる.
特徴、重み、およびラベルの特定の意味は、CERNの技術文書を参照することができる.(コンテスト公式サイトにリンクあり)
三、コード説明
1、データの読み取り
#coding:utf-8
# this is the example script to use xgboost to train
import numpy as np
import xgboost as xgb
test_size = 550000
# path to where the data lies
dpath = './data/'
# load in training data, directly use numpy delimiter skiprows converters
dtrain = np.loadtxt(dpath + '/higgsboson_training.csv', delimiter=',', skiprows=1,
converters={32: lambda x: int(x == 's'.encode('utf-8'))})
print ('finish loading from csv ')2、データの前処理
label = dtrain[:,32]
data = dtrain[:,1:31]
# rescale weight to make it same as test set
weight = dtrain[:,31] * float(test_size) / len(label) # 31 weight ，len(label) 250000
# ，
sum_wpos = sum( weight[i] for i in range(len(label)) if label[i] == 1.0 ) #
sum_wneg = sum( weight[i] for i in range(len(label)) if label[i] == 0.0 ) #
# print weight statistics
print ('weight statistics: wpos=%g, wneg=%g, ratio=%g' % ( sum_wpos, sum_wneg, sum_wneg/sum_wpos ))# 3、訓練環境の準備
トレーニングパラメータの設定:1.objective[デフォルトreg:linear]:最小化する必要がある損失関数を定義します:binary:logistic二分類の論理回帰、予測を返す確率(カテゴリではありません).Higgs Bosonコンテストは2種類の分類タスクであり,2分類の論理回帰を採用している.2. scale_pos_Weight[デフォルト1]:各種類の別サンプルが非常にアンバランスである場合、このパラメータを正の値に設定すると、アルゴリズムをより速く収束させることができる.Higgs Bosonコンテストでは,各(正/負)サンプルの重みが与えられ,すべての正/負サンプルの重みが加算され,トレーニングセットの正負サンプルの割合が得られる.3.eta[デフォルト0.3]は学習率です.オーバーフィットを防止するために、更新中に使用される収縮ステップ長.計算を上げるたびに、アルゴリズムは新しいフィーチャーの重みを直接取得します.etaは,特徴の重みを縮小することによって,リフト計算過程をより保守的にする.値の範囲は[0,1]4である.max_depth[デフォルト6]は、ツリーの最大深さを定義します.この値は、オーバーフィットを回避するためにも使用されます.max_depthが大きいほどモデルが複雑になり,より具体的で局所的なサンプルを学ぶ.標準:3-10 5.eval_metric[デフォルト値はobjectiveパラメータの値に依存する]有効なデータのメトリックメソッド.回帰問題の場合、デフォルト値はrmse、分類問題の場合、デフォルト値はerrorです.典型的な値は、rmse平均二乗誤差mae平均絶対誤差logloss負対数尤度関数値error二分類誤り率(閾値0.5)merror多分類誤り率mlogloss多分類logloss損失関数auc曲線下面積(Area Under Curve):異なるしきい値でのモデルの正解率である.6.nthread[デフォルト値は最大可能なスレッド数]このパラメータはマルチスレッド制御に使用され、システムのコア数を入力する必要があります.
CPUのすべてのコアを使用したい場合は、このパラメータを入力しないでください.アルゴリズムは自動的に検出します.
トレーニングデータはDMatixにインポートされ、後続のトレーニングがより速くなります.
# construct xgboost.DMatrix from numpy array, treat -999.0 as missing value
dtrain = xgb.DMatrix( data, label=label, missing = -999.0, weight=weight )
# setup parameters for xgboost
param = {}
# use logistic regression loss, use raw prediction before logistic transformation
# since we only need the rank
param['objective'] = 'binary:logitraw'
param['eta'] = 0.1
param['max_depth'] = 6
param['silent'] = 14、xgboostモデルトレーニング、保存
train関数のパラメータ:xgboost.train(params,dtrain,num_boost_round=10,evals=(),obj=None,feval=None,maximize=False,early_stopping_rounds=None, evals_result=None,verbose_eval=True,learning_rates=None,xgb_model=None)
params:これは辞書で、訓練中のパラメータキーワードと対応する値が含まれています.形式はparams={‘booster’:’gbtree’,’eta’:0.1}です.
dtrainトレーニングのデータ
num_boost_roundこれはリフト反復の個数を指します
evals:トレーニング中にリストの要素を評価するためのリストです.形式はevals=[(dtrain,’train’),(dval,’val’)]またはevals=[(dtrain,’train’)]である.
このコードでは、トレーニング中に検証セットの効果を観察できるようにする最初のケースを使用します.
num_round = 10
# num_round = 1000
print ('running cross validation, with preprocessing function')
# define the preprocessing function
# used to return the preprocessed training, test data, and parameter
# we can use this to do weight rescale, etc.
# as a example, we try to set scale_pos_weight
def fpreproc(dtrain, dtest, param):
label = dtrain.get_label()
ratio = float(np.sum(label == 0)) / np.sum(label==1)
param['scale_pos_weight'] = ratio
wtrain = dtrain.get_weight()
wtest = dtest.get_weight()
sum_weight = sum(wtrain) + sum(wtest)
wtrain *= sum_weight / sum(wtrain)
wtest *= sum_weight / sum(wtest)
dtrain.set_weight(wtrain)
dtest.set_weight(wtest)
return (dtrain, dtest, param)
# do cross validation, for each fold
# the dtrain, dtest, param will be passed into fpreproc
# then the return value of fpreproc will be used to generate
# results of that fold
cvresult = xgb.cv(param, dtrain, num_round, nfold=2,
metrics={'[email protected]', 'auc'}, early_stopping_rounds=2, seed = 0, fpreproc = fpreproc)
print ('finish cross validation')
print (cvresult)n_estimators = cvresult.shape[0]. # 10
from matplotlib import pyplot
# plot
test_means = cvresult['[email protected]']
test_stds = cvresult['[email protected]']
train_means = cvresult['[email protected]']
train_stds = cvresult['[email protected]']
x_axis = range(0, n_estimators)
pyplot.errorbar(x_axis, test_means, yerr=test_stds ,label='Test')
pyplot.errorbar(x_axis, train_means, yerr=train_stds ,label='Train')
pyplot.title("HiggsBoson n_estimators vs [email protected]")
pyplot.xlabel( 'n_estimators' )
pyplot.ylabel( '[email protected]' )
pyplot.savefig( 'HiggsBoson_estimators.png' ). #
#pyplot.show() # ， ， #Fit the algorithm on the data, cv refit
#alg.fit(X_train, y_train, eval_metric='[email protected]')
print ('train model using the best parameters by cv ... ')
bst = xgb.train( param, dtrain, n_estimators );
# save out model
bst.save_model('higgs_cv.model') #
print ('train finished')
pyplot.show()
5、テストプロセス
テストデータの読み取りと訓練されたモデル
#coding:utf-8
import numpy as np
import xgboost as xgb
# path to where the data lies
dpath = './data/'
modelfile = 'higgs_cv.model' #
outfile = 'higgs.pred.csv' #
# make top 15% as positive
threshold_ratio = 0.15
# load in test, directly use numpy
dtest = np.loadtxt( dpath+'/higgsboson_test.csv', delimiter=',', skiprows=1 )
data = dtest[:,1:31]
idx = dtest[:,0]
print ('finish loading from csv ')XGBoost試験環境準備
テストデータインポートDMatrixモデルインポート
xgmat = xgb.DMatrix( data, missing = -999.0 )
bst = xgb.Booster({'nthread':8}, model_file = modelfile)テスト
ypred = bst.predict( xgmat )テスト結果を整理し、結果を書き込んでファイルを提出する
res = [ ( int(idx[i]), ypred[i] ) for i in range(len(ypred)) ]
rorder = {}
for k, v in sorted( res, key = lambda x:-x[1] ):
rorder[ k ] = len(rorder) + 1
# write out predictions
ntop = int( threshold_ratio * len(rorder ) ) #threshold ，
fo = open(outfile, 'w')
nhit = 0
ntot = 0
fo.write('EventId,RankOrder,Class
') # title
for k, v in res:
if rorder[k] <= ntop:
lb = 's' #s
nhit += 1
else:
lb = 'b' #b
# change output rank order to follow Kaggle convention
fo.write('%s,%d,%s
' % ( k, len(rorder)+1-rorder[k], lb ) )
ntot += 1
fo.close()
print ('finished writing into prediction file')四、コード整理
モデル構築コード
#coding:utf-8
# this is the example script to use xgboost to train
import numpy as np
import xgboost as xgb
test_size = 550000
# path to where the data lies
dpath = './data/'
# load in training data, directly use numpy
dtrain = np.loadtxt(dpath + '/higgsboson_training.csv', delimiter=',', skiprows=1,
converters={32: lambda x: int(x == 's'.encode('utf-8'))})
print ('finish loading from csv ')
label = dtrain[:,32]
data = dtrain[:,1:31]
# rescale weight to make it same as test set
weight = dtrain[:,31] * float(test_size) / len(label)
# ，
sum_wpos = sum( weight[i] for i in range(len(label)) if label[i] == 1.0 )
sum_wneg = sum( weight[i] for i in range(len(label)) if label[i] == 0.0 )
# print weight statistics
print ('weight statistics: wpos=%g, wneg=%g, ratio=%g' % ( sum_wpos, sum_wneg, sum_wneg/sum_wpos ))
# construct xgboost.DMatrix from numpy array, treat -999.0 as missing value
dtrain = xgb.DMatrix( data, label=label, missing = -999.0, weight=weight )
# setup parameters for xgboost
param = {}
# use logistic regression loss, use raw prediction before logistic transformation
# since we only need the rank
param['objective'] = 'binary:logitraw'
param['eta'] = 0.1
param['max_depth'] = 6
param['silent'] = 1
# boost 1000 trees
num_round = 10
# num_round = 1000
print ('running cross validation, with preprocessing function')
# define the preprocessing function
# used to return the preprocessed training, test data, and parameter
# we can use this to do weight rescale, etc.
# as a example, we try to set scale_pos_weight
def fpreproc(dtrain, dtest, param):
label = dtrain.get_label()
ratio = float(np.sum(label == 0)) / np.sum(label==1)
param['scale_pos_weight'] = ratio
wtrain = dtrain.get_weight()
wtest = dtest.get_weight()
sum_weight = sum(wtrain) + sum(wtest)
wtrain *= sum_weight / sum(wtrain)
wtest *= sum_weight / sum(wtest)
dtrain.set_weight(wtrain)
dtest.set_weight(wtest)
return (dtrain, dtest, param)
# do cross validation, for each fold
# the dtrain, dtest, param will be passed into fpreproc
# then the return value of fpreproc will be used to generate
# results of that fold
cvresult = xgb.cv(param, dtrain, num_round, nfold=2,
metrics={'[email protected]', 'auc'}, early_stopping_rounds=2, seed = 0, fpreproc = fpreproc)
print ('finish cross validation')
print (cvresult)
n_estimators = cvresult.shape[0]
from matplotlib import pyplot
# plot
test_means = cvresult['[email protected]']
test_stds = cvresult['[email protected]']
train_means = cvresult['[email protected]']
train_stds = cvresult['[email protected]']
x_axis = range(0, n_estimators)
pyplot.errorbar(x_axis, test_means, yerr=test_stds, label='Test')
pyplot.errorbar(x_axis, train_means, yerr=train_stds, label='Train')
pyplot.title("HiggsBoson n_estimators vs [email protected]")
pyplot.xlabel('n_estimators')
pyplot.ylabel('[email protected]')
pyplot.savefig('HiggsBoson_estimators.png')
# pyplot.show()
#Fit the algorithm on the data, cv refit
#alg.fit(X_train, y_train, eval_metric='[email protected]')
print ('train model using the best parameters by cv ... ')
bst = xgb.train( param, dtrain, n_estimators );
# save out model
bst.save_model('higgs_cv.model')
print ('train finished')
pyplot.show()
テストコード
#coding:utf-8
import numpy as np
import xgboost as xgb
# path to where the data lies
dpath = './data/'
modelfile = 'higgs_cv.model'
outfile = 'higgs.pred.csv'
# make top 15% as positive
threshold_ratio = 0.15
# load in test, directly use numpy
dtest = np.loadtxt( dpath+'/higgsboson_test.csv', delimiter=',', skiprows=1 )
data = dtest[:,1:31]
idx = dtest[:,0]
print ('finish loading from csv ')
xgmat = xgb.DMatrix( data, missing = -999.0 )
bst = xgb.Booster({'nthread':8}, model_file = modelfile)
ypred = bst.predict( xgmat )
res = [ ( int(idx[i]), ypred[i] ) for i in range(len(ypred)) ]
rorder = {}
for k, v in sorted( res, key = lambda x:-x[1] ):
rorder[ k ] = len(rorder) + 1
# write out predictions
ntop = int( threshold_ratio * len(rorder ) )
fo = open(outfile, 'w')
nhit = 0
ntot = 0
fo.write('EventId,RankOrder,Class
')
for k, v in res:
if rorder[k] <= ntop:
lb = 's'
nhit += 1
else:
lb = 'b'
# change output rank order to follow Kaggle convention
fo.write('%s,%d,%s
' % ( k, len(rorder)+1-rorder[k], lb ) )
ntot += 1
fo.close()
print ('finished writing into prediction file')