JP6954346B2 - パラメータ推定装置、パラメータ推定方法、及びプログラム - Google Patents

パラメータ推定装置、パラメータ推定方法、及びプログラム Download PDF

Info

Publication number
JP6954346B2
JP6954346B2 JP2019515007A JP2019515007A JP6954346B2 JP 6954346 B2 JP6954346 B2 JP 6954346B2 JP 2019515007 A JP2019515007 A JP 2019515007A JP 2019515007 A JP2019515007 A JP 2019515007A JP 6954346 B2 JP6954346 B2 JP 6954346B2
Authority
JP
Japan
Prior art keywords
parameter
data points
estimation
calculated
threshold
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
JP2019515007A
Other languages
English (en)
Other versions
JPWO2018198298A1 (ja
Inventor
中野 学
学 中野
剛志 柴田
剛志 柴田
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
NEC Corp
Original Assignee
NEC Corp
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by NEC Corp filed Critical NEC Corp
Publication of JPWO2018198298A1 publication Critical patent/JPWO2018198298A1/ja
Application granted granted Critical
Publication of JP6954346B2 publication Critical patent/JP6954346B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/18Complex mathematical operations for evaluating statistical data, e.g. average values, frequency distributions, probability functions, regression analysis
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F17/00Digital computing or data processing equipment or methods, specially adapted for specific functions
    • G06F17/10Complex mathematical operations
    • G06F17/17Function evaluation by approximation methods, e.g. inter- or extrapolation, smoothing, least mean square method
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F7/00Methods or arrangements for processing data by operating upon the order or content of the data handled
    • G06F7/60Methods or arrangements for performing computations using a digital non-denominational number representation, i.e. number representation without radix; Computing devices using combinations of denominational and non-denominational quantity representations, e.g. using difunction pulse trains, STEELE computers, phase computers
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F7/00Methods or arrangements for processing data by operating upon the order or content of the data handled
    • G06F7/76Arrangements for rearranging, permuting or selecting data according to predetermined rules, independently of the content of the data

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Pure & Applied Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Mathematical Analysis (AREA)
  • Mathematical Optimization (AREA)
  • Mathematical Physics (AREA)
  • General Engineering & Computer Science (AREA)
  • Computational Mathematics (AREA)
  • Software Systems (AREA)
  • Databases & Information Systems (AREA)
  • Algebra (AREA)
  • Bioinformatics & Computational Biology (AREA)
  • Evolutionary Biology (AREA)
  • Operations Research (AREA)
  • Probability & Statistics with Applications (AREA)
  • Bioinformatics & Cheminformatics (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Computing Systems (AREA)
  • Complex Calculations (AREA)
  • Image Analysis (AREA)

Description

本発明は、外れ値を含む入力データ点に適合する幾何パラメータを推定するための、パラメータ推定装置、及びパラメータ推定方法に関し、更には、これらを実現するためのプログラムに関する。
複数のデータ点に適合する幾何パラメータを推定することは、回帰分析、又はパラメータフィッティングなどと呼ばれ、統計学、信号処理、画像処理、コンピュータビジョンといった幅広い分野における重要な要素技術となっている。推定するパラメータとしては、例えば、データ点が2次元であれば直線、放物線、ロジスティック曲線等が挙げられ、データ点が3次元であれば平面、曲面、球等が挙げられ、更に、データ点が多次元であれば超平面等が挙げられる。以下では、パラメータの種類は特に限定せず、一般的な場合について説明する。
入力となるデータ点は一般に観測ノイズを含むため、求めるべきパラメータ上の値とは完全には適合しない。これは、パラメータとデータ点との間に誤差(または「距離」とも呼ばれる)が生じているからである。そして、観測ノイズが期待値0の正規分布に従うと予想されるとき、全データ点の誤差を最小にするパラメータを推定する方法としては、最小二乗法が広く用いられている。最小二乗法が用いられる理由は、データ点の数が多ければ、観測ノイズは正規分布に近似できるという中心極限定理が成立すると考えられるからである。なお、分布が既知である限定された状況においては、データ点の事前に近似できる分布は正規分布に限らないが、以下では正規分布として説明する。
実際のパラメータ推定においては、データ点の中で、他のデータ点から大きく外れた値を持つアウトライア、(outlier、外れ値とも呼ばれる)が観測されることがある。一方で、外れ値ではないデータ点、つまり、それが含む観測ノイズが正規分布に従っており、パラメータ上の値に正しく適合する点を、インライア(inlier)と呼ぶ。アウトライアが生じる原因は、例えば、観測ミスや記録ミス、または正規分布に当てはまらない特異な個体などが原因である。最小二乗法は、全データ点が正規分布という前提のため、このようなアウトライアが1点でもあると、その大きな誤差を最小化しようと働く。また、この結果、その他のインライアに対して当てはまりの精度が低いパラメータが推定されてしまう。
このようなアウトライアを除去してパラメータを高精度に推定する方法として、非特許文献1に開示されているLMeds (Least Median Squares)、非特許文献2に開示されているRANSAC (Random Sample Consensus)、非特許文献3に開示されているM推定(M-estimator)といった推定方法が広く知られている。
非特許文献1に開示されているLMedsは、データ点からパラメータ推定に必要な最小点数をランダムに抽出して、パラメータ推定を行い、全点の誤差と誤差の中央値とを計算する。そして、LMedsは、ランダムサンプリングを一定数繰り返した後、誤差の中央値が最小となるパラメータを正しい推定とみなす。LMedsによれば、アウトライアが最大で50%含まれていても対応可能である。
また、非特許文献2に開示されているRANSACも、LMedsと同様に、データ点からパラメータ推定に必要な最小点数をランダムに抽出する。但し、RANSACは、事前に決定した閾値以下の誤差であるインライアを数え上げる。そして、RANSACは、ランダムサンプリングとインライアの数え上げとを一定数繰り返した後、インライアの数が最大となるパラメータを正しい推定とみなす。RANSACは、アウトライアが50%以上であっても対応でき、計算量がLMedsよりも少ないのが利点である。
M推定は、最小二乗法における二乗誤差の代わりに、データ点の誤差が大きいほど寄与率が小さくなるような重み付け関数を用いて、誤差の計算と重みの更新とを、誤差が収束するまで反復する。なお、非特許文献4によると、M推定はIRLS (Iterative Reweighted Least Squares:反復再重み付け最小二乗法)と数学的に等価なことが示されている。重み付けの大きさを制御するための閾値としては、M推定でも、LMeds及びRANSACのように固定値が用いられる。つまり、誤差が閾値より大きければ重みを小さくし、閾値より小さければ重みを大きくする。IRLSでは、閾値を十分大きな値で初期化し、反復のたびに0より大きく1より小さい減衰定数を掛けて徐々に閾値を減少させる。M推定は、IRLSにおいて減衰定数を1に固定した場合に等しいため、以下において、M推定の代わりにIRLSについて説明する。
P. J Rousseeuw, "Least Median of Squares Regression". Journal of the American Statistical Association. 79: 871-880, 1984. M. A. Fischler and R. C. Bolles, "Random Sample Consensus: A Paradigm for Model Fitting with Applications to Image Analysis and Automated Cartography," Comm. of the ACM, vol.24, no.6, pp.381-395, 1981. Z. Zhang, "Parameter estimation techniques: A tutorial with application to conic fitting," International journal of Image and Vision Computing, vol.15, no.1, pp.59-76, 1997. M. Black and A. Rangarajan, "On the unification of line processes, outlier rejection, and robust statistics with applications in early vision," International journal of computer vision, vol.19, no.1, pp.57-91, 1996.
しかしながら、LMeds及びRANSACといったランダムサンプリングに基づく手法には、以下のような問題がある。以下に具体的に説明する。
まず、必要な繰り返し数がアウトライアの割合に依存して指数関数的に増大する。そのため、新規のデータに対して事前に計算量を見積もることが難しい。特に、アウトライア率が50%を超えるような場合、LMedsは原理的に対応できず、RANSACは対応できるが常に高速性が要求される応用には適さない。
次に、ランダムサンプリングは原理的に同一の結果を得られず、同じデータ点であっても複数回実行した場合は結果が異なるため、安定性が低い。そして、通常はアウトライアの誤差分布が不明なため、閾値の適切な決定が困難である。閾値が大きすぎるとアウトライアをインライアと誤判定してしまう。小さすぎるとインライアをアウトライアと誤判定してしまう上に、反復回数がさらに増加する一因となる。
一方、IRLSの反復回数は、閾値の初期値と、減衰定数と、事前に定義した最小の閾値と、により決定されるため、アウトライアの割合に関係なく最大の計算量を見積もることは可能である。また、ランダムサンプリングを行わないため、、IRLSは、同一のデータ点に対しては同一の結果を返す。このように、IRLSを用いれば、上述したLMedsとRANSACとの問題の一部は解消される。
但し、IRLSを用いた場合であっても、アウトライアの誤差分布は不明であるため、閾値と減衰定数とを適切に決定することは困難である。更に、IRLSにおいては、反復毎に閾値が連続的に変更されていくため、LMedsとRANSACとの場合よりも、閾値の影響はより大きなものとなる。
本発明の目的の一例は、上記問題を解消し、複数のデータ点に適合する幾何パラメータの推定において、アウトライアとインライアとを分離する際の閾値を適切に決定し得る、パラメータ推定装置、パラメータ推定方法、およびプログラムを提供することである。
上記目的を達成するため、本発明の一側面におけるパラメータ推定装置は、
複数のデータ点をアウトライアとインライアとに分離するための閾値を計算して、前記インライアに適合するパラメータを推定するための装置であって、
前記複数のデータ点を入力として、前記パラメータを推定する、パラメータ推定部と、
前記複数のデータ点それぞれの残差の統計情報に基づいて、前記閾値を計算する、閾値決定部と、
推定された前記パラメータと計算された前記閾値とに基づいて、前記パラメータの推定が収束しているかどうかを判定し、収束していないと判定する場合は、前記パラメータ推定部、及び前記閾値決定部それぞれに、再度処理を実行させる、収束判定部と、
を備えている、ことを特徴とする。
また、上記目的を達成するため、本発明の一側面におけるパラメータ推定方法は、
複数のデータ点をアウトライアとインライアとに分離するための閾値を計算して、前記インライアに適合するパラメータを推定するための方法であって、
(a)前記複数のデータ点を入力として、前記パラメータを推定する、ステップと、
(b)前記複数のデータ点それぞれの残差の統計情報に基づいて、前記閾値を計算する、ステップと、
(c)推定された前記パラメータと計算された前記閾値とに基づいて、前記パラメータの推定が収束しているかどうかを判定する、ステップと、を有し、
前記(c)のステップで、収束していないと判定する場合は、前記(a)のステップ、及び前記(b)のステップそれぞれが再度実行される、
ことを特徴とする。
更に、上記目的を達成するため、本発明の一側面におけるプログラムは、
コンピュータによって、複数のデータ点をアウトライアとインライアとに分離するための閾値を計算して、前記インライアに適合するパラメータを推定するためのプログラムであって、
記コンピュータに、
(a)前記複数のデータ点を入力として、前記パラメータを推定する、ステップと、
(b)前記複数のデータ点それぞれの残差の統計情報に基づいて、前記閾値を計算する、ステップと、
(c)推定された前記パラメータと計算された前記閾値とに基づいて、前記パラメータの推定が収束しているかどうかを判定する、ステップと、
を実行させ、
前記(c)のステップで、収束していないと判定する場合は、前記(a)のステップ、及び前記(b)のステップそれぞれを再度実行する、
ことを特徴とする。
以上のように本発明によれば、複数のデータ点に適合する幾何パラメータの推定において、アウトライアとインライアとを分離する際の閾値を適切に決定することができる。
図1は、本発明の実施の形態のけるパラメータ推定装置の概略構成を示すブロック図である。 図2は、本発明の実施の形態におけるパラメータ推定装置の具体的構成を示すブロック図である。 図3は、本発明の実施の形態におけるパラメータ推定装置の動作を示すフロー図である。 図4は、本発明の実施の形態におけるインライアとアウトライアとの残差分布の遷移を示す図である。 図5は、本発明の実施の形態におけるパラメータ推定装置を実現するコンピュータの一例を示すブロック図である。
(実施の形態)
以下、本発明の実施の形態における、パラメータ推定装置、パラメータ推定方法、プログラムについて、図1〜図5を参照しながら説明する。
[装置構成]
最初に、図1を用いて、本実施の形態におけるパラメータ推定装置10の構成について説明する。図1は、本発明の実施の形態のけるパラメータ推定装置の概略構成を示すブロック図である。
図1に示す、本実施の形態におけるパラメータ推定装置10は、複数のデータ点をアウトライアとインライアとに分離するための閾値を計算して、インライアに適合するパラメータを推定する装置である。図1に示すように、パラメータ推定装置10は、パラメータ推定部13と、閾値決定部14と、収束判定部15とを備えている。
パラメータ推定部13は、複数のデータ点を入力として、パラメータを推定する。閾値決定部14は、複数のデータ点それぞれの残差の統計情報に基づいて、閾値を計算する。収束判定部15は、推定されたパラメータと計算された閾値とに基づいて、パラメータの推定が収束しているかどうかを判定する。また、収束判定部15は、収束していないと判定する場合は、パラメータ推定部13、及び閾値決定部14それぞれに、再度処理を実行させる。
このように、本実施の形態では、複数のデータ点をアウトライアとインライアとに分離するための閾値は、パラメータを反復して推定する際において、残差の統計情報に基づいて計算されている。このため、本実施形態によれば、複数のデータ点に適合する幾何パラメータの推定において、アウトライアとインライアとを分離する際の閾値、具体的には、IRLSにおける閾値を適切に自動的に決定することができる。
続いて、図2を用いて、本実施の形態におけるパラメータ推定装置10の構成をより具体的に説明する。図2は、本発明の実施の形態におけるパラメータ推定装置の具体的構成を示すブロック図である。図2に示すように、本実施の形態では、パラメータ推定装置10は、パラメータ推定部13、閾値決定部14、及び収束判定部15に加えて、残差計算部11と、重み更新部12とを更に備えている。
残差計算部11は、まず、入力として、複数のデータ点とパラメータの初期値とを取得する。また、残差計算部11は、取得した初期値に対する複数のデータ点それぞれの残差を計算する。
重み更新部12は、計算された残差に基づいて、複数のデータ点それぞれに設定されている重みを更新する。また、本実施の形態では、重み更新部12は、事前に決定された重み関数を用いて、各データ点の重みを更新する。
パラメータ推定部13は、本実施の形態では、複数のデータ点に加えて、更新された複数のデータ点それぞれの重みも入力として、重み付のデータ点に適合するパラメータを推定する。
閾値決定部14は、本実施の形態では、残差計算部11が計算した残差のうちインライアの残差を特定し、特定したインライアの残差の統計情報に基づいて、閾値を計算する。具体的には、閾値決定部14は、まず、閾値の初期値を取得する。次いで、閾値決定部14は、取得した初期値を、統計情報を用いて更新して、新たな閾値を計算する。その後、再度の閾値の計算が必要となると、閾値決定部14は、前回更新された閾値を更に更新する。
また、本実施の形態において、統計情報としては、インライアの残差の平均値、標準偏差、中央値、歪度、尖度、エントロピー等が挙げられる。統計情報は、インライアの残差が持つと予想される確率分布に応じて、異なる値となる。
収束判定部15は、本実施の形態では、収束していると判定する場合は、推定されたパラメータを、インライアに適合するパラメータとして外部に出力する。また、収束判定部15は、収束していないと判定する場合は、推定されたパラメータを初期値として、残差計算部11、重み更新部12、パラメータ推定部13、及び閾値決定部14それぞれに、再度処理を実行させる。つまり、本実施の形態では、パラメータの推定が収束するまで、パラメータの初期値を更新しながら、反復してパラメータの推定が行なわれる。なお、収束判定部15における判定処理の具体例については後述する。
このように、本実施の形態において、パラメータ推定装置10は、データ点をインライアとアウトライアに分離する閾値をパラメータに対するデータ点の残差に基づき決定し、更に、決定を行なう処理を、残差が収束するまで反復する。そして、収束すると、パラメータ推定装置10は、インライアに適合するパラメータを出力する。
[装置構成]
次に、本実施の形態におけるパラメータ推定装置10の動作について図3を用いて説明する。図3は、本発明の実施の形態におけるパラメータ推定装置の動作を示すフロー図である。以下の説明においては、適宜図1及び図2を参酌する。また、本実施の形態では、パラメータ推定装置10を動作させることによって、パラメータ推定方法が実施される。よって、本実施の形態におけるパラメータ推定方法の説明は、以下のパラメータ推定装置10の動作説明に代える。
図3に示すように、最初に、パラメータ推定装置10に、複数のデータ点と、パラメータの初期値と、閾値の初期とが入力されると、残差計算部11は、データ点毎に、パラメータの初期値に対する残差を計算する(ステップS1)。
ステップS1において、パラメータの初期値としては、例えば、乱数が与えられても良いし、事前に決定されているデフォルト値(例えば、直線であれば、「傾き1」、「切片0」)が与えられても良い。また、データ点とパラメータとが時系列により連続的に変化するときは、前の時刻において計算されたパラメータが初期値として入力されてもよい。また、画像特徴点のマッチングスコアように、データ点の信頼度が事前に得られている場合は、信頼度による重み付けを用いた最小二乗法によって初期値が与えられても良い。
また、閾値の初期値としては、例えば、プログラムにおいて取りうる最大の実数値が与えられても良いし、残差の上限が実験的または経験的に知られている場合は、それらが与えられていても良い。また、データ点とパラメータとが時系列により連続的に変化するときは、前の時刻における残差の最大値を定数倍して得られた値が与えられていても良い。
更に、入力されるデータ点の数は、推定すべきパラメータを表現するために必要最小限の数以上であれば良い。例えば、パラメータが直線であれば2点以上、平面であれば3点以上である。また、もし点の数がパラメータの推定に足りない場合は、パラメータ推定装置10は、実行不可能であることを示すフラグ、または実行不可能であることを示すパラメータ(例えば、全ての値がゼロまたは無限大)を出力して、処理を終了しても良い。
また、本実施の形態において、残差は、パラメータに対する各データ点の適合度合いを表す指標である。残差としては、例えば、データ点からパラメータが表す線又は面に下ろした垂線がなす幾何学的距離、データ点のノイズ分散を考慮したマハラノビス距離等が挙げられる。
次いで、重み更新部12は、ステップS1で計算された残差を入力として、事前に決定されている重み関数に従って、各データ点の重みを更新する(ステップS2)。ステップS2で用いられる重み関数としては、残差に反比例した値を返却する関数、具体的には、残差が大きいほど小さい値を返し、残差が小さいほど大きい値を返す関数が挙げられる。
次いで、パラメータ推定部13は、ステップS2で更新された各データ点の重みと、データ点とを入力として、重み付きデータ点に適合するパラメータを推定する(ステップS3)。
次いで、閾値決定部14は、ステップS1で計算された残差のうち、インライアの残差を入力として、その統計情報を利用して、新たな閾値を決定する(ステップS4)。
次に、収束判定部15は、ステップS3で推定されたパラメータと、ステップS4で計算された閾値とに基づいて、パラメータの推定が収束しているかどうかを判定する(ステップS5)。
ステップS5の判定の結果、収束していない場合は、収束判定部15は、ステップS3で推定されたパラメータをパラメータの初期値とし、ステップS4で計算された閾値を閾値の初期値とする。その後、残差計算部11が、再度、新たな初期値を用いて、再度ステップS1を実行する。これにより、ステップS2〜S4も再度実行される。
一方、ステップS5の判定の結果、収束している場合は、収束判定部15は、ステップS3で推定されたパラメータを出力する(ステップS6)。これにより、パラメータ推定装置10における処理は終了する。
また、本実施の形態における「パラメータの推定の収束」とは、一般の最適化計算又は数値計算において用いられる「収束」と同様の意味であり、これ以上反復計算を繰り返してもパラメータ及び閾値の値が変化しない、と判断される状態にあることをいう。
従って、収束判定部15は、例えば、一つ前のステップS1〜S4までの処理と最新のステップS1〜S4までの処理とにおける、パラメータの差分量及び閾値の差分量のうちの一方又は両方を計算し、差分量が一定値以下であるときに収束したと判定することができる。また、収束判定部15は、例えば、最新のステップS4で決定された閾値が一定値以下であるときに、収束したと判定することもできる。
[具体例]
続いて、パラメータ推定装置10における処理の具体例について以下に説明する。なお、以下の具体例においては、推定するパラメータをベクトルθ、i番目のデータ点をベクトルx、パラメータθに対するi番目のデータ点の残差を計算する関数をf(xi;θ)とする。また、i番目のデータ点の重みをwi、閾値をλとし、データ点は全部でN点あるとする。なお、重みは各データ点の間で相対的なものであるため、本具体例では、0から1の間の実数(0≦ wi ≦1)とする。
まず、IRLSの定義を行う。IRLSは、以下の数1で表される最適化を、閾値λを初期値から徐々に小さくしながら繰り返す操作のことである。
Figure 0006954346
ここで、関数g(wi)は、penalty functionと呼ばれ、wiが1に近づくと0、1に近づくと0となる任意の関数である。本具体例では、以下の数2が成立するとして説明する。
Figure 0006954346
ある反復(ステップS1〜S4)において、各重みwiは、上記数1をwiについて微分してゼロとした以下の数3で表される方程式を解くことで得られる。
Figure 0006954346
これを上記数2の場合について解くと、以下の数4が得られる。
Figure 0006954346
上記数4をλの関数と見れば、0から1の値をとることが確認できる。また、上記数4はGeman-McClureによるM推定の重み付け関数と一致することが、上述した非特許文献4に示されている。
ここで、上述した具体例を、図3に示したステップS1〜S5に沿って更に説明する。まず、ステップS1において、残差計算部11は、パラメータθに対する各データ点xiの残差f(xi;θ)を計算する。
次に、ステップS2において、重み更新部12は、閾値λと、上記数4に従って、各データ点xiの重みwiを更新する。
次に、ステップS3において、パラメータ推定部13は、ステップS2で更新された重みwiを上記数1に代入して、最適化問題を最小化するパラメータθを推定する。最適化問題の解は、例えば、上記数1がmf(xi;θ)の二乗和であることを利用して、ガウスニュートン法、レーベンバーグマーカード法などによって求めることができる。
次に、ステップS4において、閾値決定部14は、ステップS1で計算された残差f(xi;θ)と、ステップS2で更新された重みwiとを用いて、新たな閾値λを決定する。
具体的には、閾値決定部14は、まず、重みwiから各データ点がインライアかどうかを判定する。重みは0から1の間の値を取るから、閾値決定部14は、例えば重みwiが0.5以上のデータ点をインライアとして判定する。なお、重みwiが1に近い値となるときにインライアとして判定するとすれば、より厳しい判定結果となる。
次に、インライアの残差だけを用いて、閾値決定部14は、新たなしきい値λを決定する。ここで、インライアの残差が、期待値μ、標準偏差σの正規分布に従うと仮定すると、インライアの残差は、ある定数βを用いてμ+βσの範囲内に分布すると考えられる。例えば、β=1、β=2とすると、残差f(xi;θ)は、それぞれ約84%、約97%となり、これらはμ+βσ以下となる。つまり、λ=μ+βσと更新すれば、その反復時点の任意の割合のインライアを抽出できる。
次に、ステップS5において、収束判定部15は、ステップS3で推定されたパラメータθとステップS4で決定された閾値λとを入力として、パラメータ推定が収束したかどうかを判定する。
上述したように、例えば、収束判定部15、前回のステップS1〜S4の実行によって得られたパラメータθ又は閾値λと、前回のステップS1〜S4の実行によって得られたパラメータθ又は閾値λとの差分を比較し、差分が一定値以下であれば収束したと判定する。また、収束判定部15は、最新のステップS4で決定された閾値λが事前に決めていた最小値よりも小さくなった場合に、収束したと判定しても良い。
判定の結果、収束していない場合は、収束判定部15は、ステップS3で推定されたパラメータをパラメータの初期値とし、ステップS4で計算された閾値を閾値の初期値とする。その後、再度ステップS1以降が実行される。一方、収束している場合は、収束判定部15は、ステップS3で推定されたパラメータを出力する
[実施の形態における効果]
以上のように、本実施の形態によれば、IRLSにおける閾値を自動的に且つ適切に決定でき、アウトライアとインライアとを高精度に分離できる。また、この結果、インライアに適合するパラメータを高精度に推定可能である。
ここで、上述したステップS1〜S4を反復して実行することにより、インライアとアウトライアとの分離ができる原理を、図4を用いて説明する。図4は、本発明の実施の形態におけるインライアとアウトライアとの残差分布の遷移を示す図である。図4は、上段から下段へと初期状態から収束状態までを示している。
図4に示すように、上段の初期状態、つまり、インライアとアウトライアの区別がほとんどつかない状態では、2つの残差分布が混じっている。中心極限定理により、これは混合正規分布と仮定できる。ここで、ステップS4においてλ=μ+βσとすることは、μ+βσより大きい残差f(xi;θ)を持つデータ点の重みwiが急激に小さくなることを意味する。そして、データ点により大きな重みが与えられることは、パラメータにおいてインライアの占める割合が上昇することと等価である。インライア率が高まった状態で、上記の数1を解くと、よりインライアに適合するパラメータθが求まり、図4の中段の状態となる。そして、ステップS1〜S4が反復されると、最終的に図4の下段に示すように、インライアとアウトライアとが分離された状態に収束する。
[変形例]
ここで、本実施の形態における変形例1〜4について説明する。本発明は、上述した実施の形態のみに限定されることはない。本発明は、上述した実施の形態に対して、いわゆる当業者が理解し得る多様な変更を適用することが可能である。例えば、本発明は、以下の変形例に示す形態によっても実施可能である。
(1)変形例1
どのpenalty functionを用いるかは、上述した実施の形態に限定されない。例えば、非特許文献4にはpenalty functionの様々な具体例が詳しく記載されている。また、非特許文献4にはIRLSとM推定との間の変換方法が記載されているため、非特許文献3に記載されている様々なM推定の関数が利用できる。
(2)変形例2
閾値決定方法は、上述した実施の形態に限定されない。例えば、前回の反復(ステップS1〜S4の実施)で得られた閾値を保存しておき、この閾値と今回得られた閾値との差分が一定値以下となるようにして、急激な減衰が抑制されても良い。また、例えば、従来の減衰定数を用いる方法も同時に計算しておき、より小さな値又はより大きな値が、新たな閾値として用いられても良い。また、例えば、インライアの残差分布がコーシー分布である場合、期待値と標準偏差とは存在しないため、代わりに最頻値と半値半幅とを与えるスケール定数が用いられても良い。
(3)変形例3
収束判定方法は、上述した実施の形態に限定されない。例えば、収束の判定基準として、反復の前後における上記数1で表される目的関数の差分が用いられても良いし、データ点の重みwiの差分が用いられても良い。このとき、収束判定部15は、上記数1で表される目的関数または重みwiを入力として受け付ける。
(4)変形例4
パラメータ推定装置10に入力されるデータは、上述の実施の形態に示された例に限定されない。例えば、パラメータの初期値の代わりに、重みの初期値が入力されても良い。この場合、パラメータ推定装置10においては、パラメータ推定部13が最初に処理を実行することになり、次いで残差計算部11、重み更新部12、閾値決定部14、収束判定部15の順で処理が実行される。なお、各部による処理の実行順序は、適宜変更可能である。
[プログラム]
本実施の形態におけるプログラムは、コンピュータに、図3に示すステップS1〜S6を実行させるプログラムであれば良い。このプログラムをコンピュータにインストールし、実行することによって、本実施の形態におけるパラメータ装置とパラメータ推定方法とを実現することができる。この場合、コンピュータのCPU(Central Processing Unit)は、残差計算部11、重み更新部12、パラメータ推定部13、閾値決定部14、及び収束判定部15として機能し、処理を行なう。
また、本実施の形態におけるプログラムは、複数のコンピュータによって構築されたコンピュータシステムによって実行されても良い。この場合は、例えば、各コンピュータが、それぞれ、残差計算部11、重み更新部12、パラメータ推定部13、閾値決定部14、及び収束判定部15のいずれかとして機能しても良い。
ここで、本実施の形態におけるプログラムを実行することによって、パラメータ推定装置10を実現するコンピュータについて図5を用いて説明する。図5は、本発明の実施の形態におけるパラメータ推定装置を実現するコンピュータの一例を示すブロック図である。
図5に示すように、コンピュータ110は、CPU111と、メインメモリ112と、記憶装置113と、入力インターフェイス114と、表示コントローラ115と、データリーダ/ライタ116と、通信インターフェイス117とを備える。これらの各部は、バス121を介して、互いにデータ通信可能に接続される。
CPU111は、記憶装置113に格納された、本実施の形態におけるプログラム(コード)をメインメモリ112に展開し、これらを所定順序で実行することにより、各種の演算を実施する。メインメモリ112は、典型的には、DRAM(Dynamic Random Access Memory)等の揮発性の記憶装置である。また、本実施の形態におけるプログラムは、コンピュータ読み取り可能な記録媒体120に格納された状態で提供される。なお、本実施の形態におけるプログラムは、通信インターフェイス117を介して接続されたインターネット上で流通するものであっても良い。
また、記憶装置113の具体例としては、ハードディスクドライブの他、フラッシュメモリ等の半導体記憶装置が挙げられる。入力インターフェイス114は、CPU111と、キーボード及びマウスといった入力機器118との間のデータ伝送を仲介する。表示コントローラ115は、ディスプレイ装置119と接続され、ディスプレイ装置119での表示を制御する。
データリーダ/ライタ116は、CPU111と記録媒体120との間のデータ伝送を仲介し、記録媒体120からのプログラムの読み出し、及びコンピュータ110における処理結果の記録媒体120への書き込みを実行する。通信インターフェイス117は、CPU111と、他のコンピュータとの間のデータ伝送を仲介する。
また、記録媒体120の具体例としては、CF(Compact Flash(登録商標))及びSD(Secure Digital)等の汎用的な半導体記憶デバイス、フレキシブルディスク(Flexible Disk)等の磁気記録媒体、又はCD−ROM(Compact Disk Read Only Memory)などの光学記録媒体が挙げられる。
なお、本実施の形態におけるパラメータ推定装置10は、プログラムがインストールされたコンピュータではなく、各部に対応したハードウェアを用いることによっても実現可能である。更に、パラメータ推定装置10は、一部がプログラムで実現され、残りの部分がハードウェアで実現されていてもよい。
上述した実施の形態の一部又は全部は、以下に記載する(付記1)〜(付記15)によって表現することができるが、以下の記載に限定されるものではない。
(付記1)
複数のデータ点をアウトライアとインライアとに分離するための閾値を計算して、前記インライアに適合するパラメータを推定するための装置であって、
前記複数のデータ点を入力として、前記パラメータを推定する、パラメータ推定部と、
前記複数のデータ点それぞれの残差の統計情報に基づいて、前記閾値を計算する、閾値決定部と、
推定された前記パラメータと計算された前記閾値とに基づいて、前記パラメータの推定が収束しているかどうかを判定し、収束していないと判定する場合は、前記パラメータ推定部、及び前記閾値決定部それぞれに、再度処理を実行させる、収束判定部と、
を備えている、ことを特徴とするパラメータ推定装置。
(付記2)
前記複数のデータ点と前記パラメータの初期値とを取得し、取得した前記初期値に対する前記複数のデータ点それぞれの残差を計算する、残差計算部と、
計算された前記残差に基づいて、前記複数のデータ点それぞれに設定されている重みを更新する、重み更新部と、を更に備え、
前記パラメータ推定部は、前記複数のデータ点に加えて、更新された前記複数のデータ点それぞれの重みも入力として、前記パラメータを推定し、
前記閾値決定部は、計算された前記残差のうち前記インライアの残差を特定し、特定した前記インライアの残差の前記統計情報に基づいて、前記閾値を計算し、
前記収束判定部は、
収束していると判定する場合は、推定された前記パラメータを、前記インライアに適合するパラメータとして出力し、
収束していないと判定する場合は、推定された前記パラメータを前記初期値として、前記残差計算部、前記重み更新部、前記パラメータ推定部、及び前記閾値決定部それぞれに、再度処理を実行させる、
付記1に記載のパラメータ推定装置。
(付記3)
前記パラメータ推定部が、反復再重み付け最小二乗法、又はM推定を実行して、前記パラメータを推定する、
付記2に記載のパラメータ推定装置。
(付記4)
前記閾値決定部は、前記インライアの残差の事前確率分布から前記統計情報を求め、求めた前記統計情報に基づいて、前記閾値を計算する、
付記2または3に記載のパラメータ推定装置。
(付記5)
前記事前確率分布は、混合正規分布で表され、
前記閾値決定部は、前記統計情報として、前記インライアの残差の平均値と標準偏差とを求める、
付記4に記載のパラメータ推定装置。
(付記6)
複数のデータ点をアウトライアとインライアとに分離するための閾値を計算して、前記インライアに適合するパラメータを推定するための方法であって、
(a)前記複数のデータ点を入力として、前記パラメータを推定する、ステップと、
(b)前記複数のデータ点それぞれの残差の統計情報に基づいて、前記閾値を計算する、ステップと、
(c)推定された前記パラメータと計算された前記閾値とに基づいて、前記パラメータの推定が収束しているかどうかを判定する、ステップと、を有し、
前記(c)のステップで、収束していないと判定する場合は、前記(a)のステップ、及び前記(b)のステップそれぞれが再度実行される、
ことを特徴とするパラメータ推定方法。
(付記7)
(d)前記複数のデータ点と前記パラメータの初期値とを取得し、取得した前記初期値に対する前記複数のデータ点それぞれの残差を計算する、ステップと、
(e)計算された前記残差に基づいて、前記複数のデータ点それぞれに設定されている重みを更新する、ステップと、を更に有し、
前記(a)のステップで、前記複数のデータ点に加えて、前記(e)のステップで更新された前記複数のデータ点それぞれの重みも入力として、前記パラメータを推定し、
前記(b)のステップで、前記(d)のステップで計算された前記残差のうち前記インライアの残差を特定し、特定した前記インライアの残差の前記統計情報に基づいて、前記閾値を計算し、
前記(c)のステップで、収束していると判定する場合は、前記(c)のステップにおいて、推定された前記パラメータを、前記インライアに適合するパラメータとして出力し、
前記(c)のステップで、収束していないと判定する場合は、前記(a)のステップで推定された前記パラメータを前記初期値として、前記(d)のステップ、前記(e)のステップ、前記(a)のステップ、及び前記(b)のステップそれぞれが再度実行される、
付記6に記載のパラメータ推定方法。
(付記8)
前記(a)のステップにおいて、反復再重み付け最小二乗法、又はM推定を実行して、前記パラメータを推定する、
付記7に記載のパラメータ推定方法。
(付記9)
前記(b)のステップにおいて、前記インライアの残差の事前確率分布から前記統計情報を求め、求めた前記統計情報に基づいて、前記閾値を計算する、
付記7または8に記載のパラメータ推定方法。
(付記10)
前記事前確率分布は、混合正規分布で表され、
前記(b)のステップにおいて、前記統計情報として、前記インライアの残差の平均値と標準偏差とを求める、
付記9に記載のパラメータ推定方法。
(付記11)
コンピュータによって、複数のデータ点をアウトライアとインライアとに分離するための閾値を計算して、前記インライアに適合するパラメータを推定するためのプログラムであって、
記コンピュータに、
(a)前記複数のデータ点を入力として、前記パラメータを推定する、ステップと、
(b)前記複数のデータ点それぞれの残差の統計情報に基づいて、前記閾値を計算する、ステップと、
(c)推定された前記パラメータと計算された前記閾値とに基づいて、前記パラメータの推定が収束しているかどうかを判定する、ステップと、
を実行させ、
前記(c)のステップで、収束していないと判定する場合は、前記(a)のステップ、及び前記(b)のステップそれぞれを再度実行する、
ことを特徴とするプログラム
(付記12)
記コンピュータに、
(d)前記複数のデータ点と前記パラメータの初期値とを取得し、取得した前記初期値に対する前記複数のデータ点それぞれの残差を計算する、ステップと、
(e)計算された前記残差に基づいて、前記複数のデータ点それぞれに設定されている重みを更新する、ステップと、
更に実行させ、
前記(a)のステップで、前記複数のデータ点に加えて、前記(e)のステップで更新された前記複数のデータ点それぞれの重みも入力として、前記パラメータを推定し、
前記(b)のステップで、前記(d)のステップで計算された前記残差のうち前記インライアの残差を特定し、特定した前記インライアの残差の前記統計情報に基づいて、前記閾値を計算し、
前記(c)のステップで、収束していると判定する場合は、前記(c)のステップにおいて、推定された前記パラメータを、前記インライアに適合するパラメータとして出力し、
前記(c)のステップで、収束していないと判定する場合は、前記(a)のステップで推定された前記パラメータを前記初期値として、前記(d)のステップ、前記(e)のステップ、前記(a)のステップ、及び前記(b)のステップそれぞれを再度実行する、
付記11に記載のプログラム
(付記13)
前記(a)のステップにおいて、反復再重み付け最小二乗法、又はM推定を実行して、前記パラメータを推定する、
付記12に記載のプログラム
(付記14)
前記(b)のステップにおいて、前記インライアの残差の事前確率分布から前記統計情報を求め、求めた前記統計情報に基づいて、前記閾値を計算する、
付記12または13に記載のプログラム
(付記15)
前記事前確率分布は、混合正規分布で表され、
前記(b)のステップにおいて、前記統計情報として、前記インライアの残差の平均値と標準偏差とを求める、
付記14に記載のプログラム
以上、実施の形態を参照して本願発明を説明したが、本願発明は上記実施の形態に限定されるものではない。本願発明の構成や詳細には、本願発明のスコープ内で当業者が理解し得る様々な変更をすることができる。
以上のように本発明によれば、複数のデータ点に適合する幾何パラメータの推定において、アウトライアとインライアとを分離する際の閾値を適切に決定することができる。本発明は、複数のデータ点に適合する幾何パラメータを推定することが求められる分野、例えば、統計学、信号処理、画像処理、コンピュータビジョン等に有用である。
10 カメラパラメータ推定装置
11 残差計算部
12 重み更新部
13 パラメータ推定部
14 閾値決定部
15 収束判定部
110 コンピュータ
111 CPU
112 メインメモリ
113 記憶装置
114 入力インターフェイス
115 表示コントローラ
116 データリーダ/ライタ
117 通信インターフェイス
118 入力機器
119 ディスプレイ装置
120 記録媒体
121 バス

Claims (5)

  1. 複数のデータ点をアウトライアとインライアとに分離するための閾値を計算して、前記インライアに適合するパラメータを推定するための装置であって、
    前記複数のデータ点と前記パラメータの初期値とを取得し、取得した前記初期値に対する前記複数のデータ点それぞれの残差を計算する、残差計算部と、
    計算された前記残差に基づいて、前記複数のデータ点それぞれに設定されている重みを更新する、重み更新部と、
    前記複数のデータ点に加えて、更新された前記複数のデータ点それぞれの重みも入力として、前記パラメータを推定する、パラメータ推定部と、
    計算された前記残差のうちインライアの残差は、期待値と標準偏差の正規分布に従うと仮定し、前記期待値と前記標準偏差と定数に基づいて前記閾値を計算する、閾値決定部と、
    推定された前記パラメータと計算された前記閾値とに基づいて、前記パラメータの推定が収束しているかどうかを判定し、収束していないと判定する場合は、推定された前記パラメータを前記初期値として、前記残差計算部、前記重み更新部、前記パラメータ推定部、及び前記閾値決定部それぞれに、再度処理を実行させる、収束判定部と、
    を備えている、ことを特徴とするパラメータ推定装置。
  2. 記収束判定部は、収束していると判定する場合は、推定された前記パラメータを、前記インライアに適合するパラメータとして出力する、
    請求項1に記載のパラメータ推定装置。
  3. 前記パラメータ推定部が、反復再重み付け最小二乗法、又はM推定を実行して、前記パラメータを推定する、
    請求項2に記載のパラメータ推定装置。
  4. 複数のデータ点をアウトライアとインライアとに分離するための閾値を計算して、前記インライアに適合するパラメータを推定するための方法であって、
    前記複数のデータ点と前記パラメータの初期値とを取得し、取得した前記初期値に対する前記複数のデータ点それぞれの残差を計算し、
    計算された前記残差に基づいて、前記複数のデータ点それぞれに設定されている重みを更新し、
    前記複数のデータ点に加えて、更新された前記複数のデータ点それぞれの重みも入力として、前記パラメータを推定し、
    計算された前記残差のうちインライアの残差は、期待値と標準偏差の正規分布に従うと仮定し、前記期待値と前記標準偏差と定数に基づいて前記閾値を計算し、
    推定された前記パラメータと計算された前記閾値とに基づいて、前記パラメータの推定が収束しているかどうかを判定し、収束していないと判定する場合は、推定された前記パラメータを前記初期値として、前記残差の計算、前記重みの更新、前記パラメータの推定、及び前記閾値の決定それぞれを、再度実行させる、
    ことを特徴とするパラメータ推定方法。
  5. コンピュータによって、複数のデータ点をアウトライアとインライアとに分離するための閾値を計算して、前記インライアに適合するパラメータを推定するためのプログラムであって、
    前記コンピュータに、
    前記複数のデータ点と前記パラメータの初期値とを取得し、取得した前記初期値に対する前記複数のデータ点それぞれの残差を計算させ、
    計算された前記残差に基づいて、前記複数のデータ点それぞれに設定されている重みを更新させ、
    前記複数のデータ点に加えて、更新された前記複数のデータ点それぞれの重みも入力として、前記パラメータを推定させ、
    計算された前記残差のうちインライアの残差は、期待値と標準偏差の正規分布に従うと仮定し、前記期待値と前記標準偏差と定数に基づいて前記閾値を計算させ、
    推定された前記パラメータと計算された前記閾値とに基づいて、前記パラメータの推定が収束しているかどうかを判定し、収束していないと判定する場合は、推定された前記パラメータを前記初期値として、前記残差の計算、前記重みの更新、前記パラメータの推定、及び前記閾値の決定それぞれを、再度実行させる、
    ことを特徴とするプログラム。
JP2019515007A 2017-04-27 2017-04-27 パラメータ推定装置、パラメータ推定方法、及びプログラム Active JP6954346B2 (ja)

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
PCT/JP2017/016869 WO2018198298A1 (ja) 2017-04-27 2017-04-27 パラメータ推定装置、パラメータ推定方法、及びコンピュータ読み取り可能な記録媒体

Publications (2)

Publication Number Publication Date
JPWO2018198298A1 JPWO2018198298A1 (ja) 2020-03-05
JP6954346B2 true JP6954346B2 (ja) 2021-10-27

Family

ID=63920210

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2019515007A Active JP6954346B2 (ja) 2017-04-27 2017-04-27 パラメータ推定装置、パラメータ推定方法、及びプログラム

Country Status (3)

Country Link
US (1) US11455372B2 (ja)
JP (1) JP6954346B2 (ja)
WO (1) WO2018198298A1 (ja)

Families Citing this family (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113658156A (zh) * 2021-08-24 2021-11-16 凌云光技术股份有限公司 一种去除深度图像中局外点的球体拟合方法和装置

Family Cites Families (4)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2005201869A (ja) 2004-01-19 2005-07-28 Mitsutoyo Corp 信号処理方法、信号処理プログラム、この信号処理プログラムを記録した記録媒体および信号処理装置
US20070255512A1 (en) * 2006-04-28 2007-11-01 Delenstarr Glenda C Methods and systems for facilitating analysis of feature extraction outputs
US8611657B2 (en) * 2011-02-25 2013-12-17 Adobe Systems Incorporated Robust fitting of surfaces from noisy data
WO2013087884A1 (en) * 2011-12-16 2013-06-20 Universiteit Gent Mutant microorganisms to synthesize colanic acid, mannosylated and/or fucosylated oligosaccharides

Also Published As

Publication number Publication date
US11455372B2 (en) 2022-09-27
JPWO2018198298A1 (ja) 2020-03-05
WO2018198298A1 (ja) 2018-11-01
US20200097522A1 (en) 2020-03-26

Similar Documents

Publication Publication Date Title
Kurtz et al. Cross-entropy-based adaptive importance sampling using Gaussian mixture
CN104869126B (zh) 一种网络入侵异常检测方法
US11030246B2 (en) Fast and accurate graphlet estimation
KR101852527B1 (ko) 기계학습 기반의 동적 시뮬레이션 파라미터 교정 방법
JP7091872B2 (ja) 検知装置及び検知方法
JP2011113550A (ja) 情報処理装置、情報処理方法、情報処理システム、プログラム及びデータ構造
CN111652371A (zh) 一种离线强化学习网络训练方法、装置、系统及存储介质
CN111078639A (zh) 数据标准化方法、装置以及电子设备
CN114424204A (zh) 使用强化学习的数据评估
US20220058312A1 (en) Estimation apparatus, optimization apparatus, estimation method, optimization method, and program
JP6954346B2 (ja) パラメータ推定装置、パラメータ推定方法、及びプログラム
CN111612022A (zh) 用于分析数据的方法、设备和计算机存储介质
US12002549B2 (en) Knowledge reuse-based method and system for predicting cell concentration in fermentation process
JP6398991B2 (ja) モデル推定装置、方法およびプログラム
CN116543259A (zh) 一种深度分类网络噪声标签建模与纠正方法、系统及存储介质
CN115272152A (zh) 一种对抗医学图像的生成方法、装置、设备及存储介质
JP6233432B2 (ja) 混合モデルの選択方法及び装置
JP7477859B2 (ja) 計算機、計算方法及びプログラム
JP7331938B2 (ja) 学習装置、推定装置、学習方法及び学習プログラム
JP7056804B2 (ja) 経験損失推定システム、経験損失推定方法および経験損失推定プログラム
CN113656669A (zh) 标签更新方法及装置
CN114693052A (zh) 风险预测模型的训练方法、装置、计算设备和介质
JP5970579B2 (ja) 混合モデル決定用の装置、方法、およびプログラム
EP3940601A1 (en) Information processing apparatus, information processing method, and information program
JP6453618B2 (ja) 算出装置、方法及びプログラム

Legal Events

Date Code Title Description
A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20191021

A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20191021

A131 Notification of reasons for refusal

Free format text: JAPANESE INTERMEDIATE CODE: A131

Effective date: 20210105

A521 Request for written amendment filed

Free format text: JAPANESE INTERMEDIATE CODE: A523

Effective date: 20210303

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20210831

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20210913

R150 Certificate of patent or registration of utility model

Ref document number: 6954346

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R150