JP6029683B2 - データ解析装置、データ解析プログラム - Google Patents

Description

本発明は、サンプルデータをクラスタリング解析する技術に関するものである。

再生医療、遺伝子を用いた診断、あるいは生命現象の基本的な理解のため、組織の平均としての遺伝子の発現量ばかりでなく、組織を構成する１つ１つの細胞の中味を定量的に分析することが重要視され始めている。このような個々の細胞の特性を１つずつを解析することを単一細胞解析と呼んでいる。

単一細胞中の生体分子の量が極めて微量であるため、単一細胞解析は、これまで細胞膜上の蛋白質など一部の生体分子を対象とする解析にしか用いられてこなかった。しかし近年では、技術の発展により単一細胞中の微量な生体分子を定量評価できるようになりつつある。

下記非特許文献１には、単一細胞中の遺伝子発現解析に関して、定量ＰＣＲ装置を用いた方法によって、十分な精度で特定の遺伝子発現量を計測できる方法が開示されている。下記非特許文献２には、同様に単一細胞中に遺伝子発現解析に関して、大規模ＤＮＡシーケンサ（次世代シーケンサ）を用いて、ほぼすべての遺伝子の発現を定量化する方法が開示されている。非特許文献２には、細胞の種類を特定するためのデータ解析方法も開示されている。今後、ゲノム配列、細胞内の蛋白質、さらに細胞内の様々な生体分子が単一細胞レベルで同定されるようになることが予想される。

Ｎａｔｕｒｅ Ｍｅｔｈｏｄ，Ｖｏｌ．６，Ｎｏ．７（２００９），ｐｐ．５０３ Ｇｅｎｏｍｅ Ｒｅｓｅａｒｃｈ，Ｖｏｌ．２１，Ｎｏ．７（２０１１），ｐｐ．１０８８

上記のような単一細胞解析技術の進歩によって、これまで均一であると仮定され解析されてきた細胞組織は、単一細胞解析で得られるデータを用いて、これまで知られているよりも詳細なグループ、下位組織を構成していることが明らかとなるであろう。その結果、非常に多数の細胞から構成されている人体などの個体中の複雑な生命現象は、細胞のデータによって分類された細胞のグループで構成され、このグループ間を様々な生化学的なシグナルがやりとりされるネットワークとして生命が把握され、生命科学分野、特に医療あるいは創薬分野において大きなインパクトを与えると考えられる。

たとえば、これまで均一であると考えられてきたがん組織をグループ分けし、各グループに対して遺伝子変異を解析することによって、より適切な分子診断薬を選択することができるようになる可能性がある。また、血中の免疫細胞中の遺伝子発現量を解析することにより、様々な病気を診断することができる可能性が示唆されているが、免疫細胞の詳細な分類によって、より精度の高い診断ができるようになると考えられる。

しかし、データのみを用いて細胞を分類するアルゴリズムおよびこれを応用した解析／診断装置は、細胞を分類して医療向けの診断において用いるために必要な特性が必ずしも十分ではない。ここでいう必要な特性の例としては、例えば細胞のグループ化（以下ではクラスタリングと呼ぶ）において、最適なグループ数（以下ではクラスタ数と呼ぶ）が事前に分かっていない場合であっても、適切なクラスタリングを実施することなどが考えられる。特に、データ数の少ない例外的なクラスタを独立のクラスタとするのか、他のデータ数の多いクラスタの一部とするのかを判定することは、従来の解析／診断装置においては困難である。

本発明は、上記のような課題に鑑みてなされたものであり、実験誤差などに起因する例外データが生じる場合においても、クラスタリングを適切に実施することができるデータ解析装置を提供することを目的とする。

本発明に係るデータ解析装置は、実験誤差データが記述している実験誤差の範囲に応じて、クラスタ境界を緩和するためのクラスタ範囲パラメータをあらかじめ定めておく。クラスタリングの過程において、いずれのクラスタにも属さない例外データについては、例外データからさらにクラスタ範囲パラメータによって定められる距離だけ離れた領域がいずれかのクラスタに含まれる場合は、例外データがそのクラスタに属するものと判定し、前記距離だけ離れてなおいずれのクラスタにも含まれない場合は、独立したクラスタを構成するものと判定する。

本発明に係るデータ解析装置によれば、単一細胞の解析結果を用いて細胞を適切に分類することができる。また、分類された細胞の種類数を精度よく決定することができる。

模擬サンプルデータのクラスタリング結果を示す図である。 階層的クラスタリング法について説明する図である。 赤池情報量基準を図１と同じ模擬サンプルデータに対して適用した結果を示す図である。 実施形態１に係るデータ解析装置がサンプルデータをクラスタリングする処理の概略フローを示す図である。 ステップＳ４０５の詳細を説明する処理フロー図である。 ステップＳ５０５〜Ｓ５０７の処理イメージを示す図である。 ステップＳ５０８〜Ｓ５０９の処理イメージを示す図である。 図１で説明した模擬サンプルデータについてＣＲ値を１、２、４に掃引し、クラスタ数の尤もらしさを表す対数尤度を計算した結果を示す。 実施形態２に係るデータ解析装置９０６の機能ブロック図である。 実施形態２において、クラスタ数の尤もらしさを表す対数尤度を計算した結果を示す。 階層的クラスタリングの結果と本実施形態２に係るデータ解析装置９０６によるクラスタリング結果を比較する図である。 各クラスタ（１〜５の番号で表示）のＢｇｌａｐ１、ＰＰａｒｇ, Ｃｏｌ２ａ１の遺伝子発現量を示す図である。 図１２で説明した３つの遺伝子について、クラスタリングを実施せず、分化誘導ありの場合となしの場合について発現量を示した図である。 実施形態３に係るデータ解析装置９０６の機能ブロック図である。 実施形態４に係るデータ解析装置９０６の機能ブロック図である。

＜従来のデータ解析手法＞
以下では本発明の理解を促進するため、まず初めに従来のデータ解析手法およびその課題について、具体例を挙げて説明する。その後に、本発明に係るデータ解析装置の具体的な構成について説明する。

単一細胞解析で得られた細胞データは、たとえば主因子分析法を用いて解析することができる。主因子分析法によって得られたデータは、目視によってグループを判定するために用いられることが多い。従来のクラスタリング方法の課題を具体的に見るため、以下では模擬サンプルデータを用いることにする。具体的には、４つの遺伝子、１８０個の細胞について単一細胞解析を実施したことを想定し、コンピュータ上で乱数を用いて細胞データを生成した。

図１は、模擬サンプルデータのクラスタリング結果を示す図である。図１（ａ）は模擬データ内の各クラスタの構成を示す。εは乱数であり、平均値０標準偏差σの分布に従うものとした。クラス２とクラスタ３の間の距離を制御するためのパラメータαを設定し、α＝６×σに設定した。乱数εは４次元空間において等方的なガウス分布に従うと仮定しており、その片側標準偏差を上記σとした。

図１（ｂ）は、図１（ａ）に示す模擬データについて主因子分析を実施し、上位２個の主因子について可視化したデータである。第１主因子に対応する軸をＰＣ１とし、第２主因子に対応する軸をＰＣ２とした。各データ集合に添えている数字は、図１（ａ）に記載している各クラスタ番号である。図１（ｂ）に示すように、６つのクラスタを容易に目視確認することができる。しかし、主因子分析法はクラスタリングを実施するためのアルゴリズムではないため、目視確認以外の手段によって明確にクラスタリング結果を与えるものではない。また、クラスタリング結果の信頼度も定量的に与えられるものではない。

クラスタリング結果の信頼度は、複数の被験者から得た細胞データを比較する際に用いることができる。例えば、健常な被験者から採取した細胞データのクラスタリング結果が信頼できる場合、そのクラスタリング結果を基準として他の被験者の細胞データを検証することができる。しかしクラスタリング結果の信頼度が不明であれば、これを基準として用いることが適切であるのか否か判断することができない。したがって、細胞解析においてクラスタリング結果の信頼度を得ることは、非常に有用である。

図２は、階層的クラスタリング法について説明する図である。階層的クラスタリング法は、主因子分析法と同様にデータを分類する手法としてよく用いられる。ここでは、図１と同じデータを階層的クラスタリング法によってクラスタリングした結果を示す。各細胞のデータは４次元ユークリッド空間内のベクトルとして表現されているため、各細胞間の距離はユークリッド距離を用いて評価した。

階層的クラスタリング法においては、データ間の距離が最も近い２つのデータをペアにして、その距離に応じた高さのトーナメント図に現れる矩形を設定する。以下、ペアの平均値の位置に２つのデータを代表するデータがあるものと仮定して次のデータに進み、すべてのデータがペアとして結びつけられるまで同様の処理を繰り返す。トーナメント図における矩形の高さが高いほど、クラスタ間の距離が大きいことを示す。ある高さにおいてトーナメント図を水平方向に切断した時に得られる長い鉛直線との交点の数がクラスタ数に相当すると考えられる。

しかし、階層的クラスタリング法においては、最適なクラスタ数は求めることができない。図２に示す例においては、５個のクラスタあるいは６個のクラスタいずれも妥当なクラスタ数であるように思われる。すなわち、階層的クラスタリング法によって得られるクラスタ数は、主因子分析法において目視確認により得られるクラスタ数とは必ずしも一致しない。したがって、階層的クラスタリング法は、最適なクラスタ数およびその時のクラスタリングの信頼度を求めることは困難である。

クラスタリングの信頼度を評価するためには、サンプルデータを確率分布へフィッティングさせることが最も自然である。確率分布へフィッティングさせる手法としては、混合ガウシアンモデルと呼ばれる方法が最もよく知られている。混合ガウシアンモデルにおいては、各細胞データはガウス分布していると仮定し、各細胞データはいずれかのクラスタに所属していると見なす。次に、ガウシアン確率密度関数に対する対数尤度を計算し、最尤法によってクラスタ数、クラスタ位置（平均値）、分布（標準偏差）を決定する。

一般に、単純に対数尤度を用いてクラスタ数などを決定する場合においては、クラスタ数とデータ数が一致するまでクラスタ総数を増やせば良いフィッティングが得られる。したがって、単一細胞解析のように最適なクラスタ総数を決定したい場合においては不向きである。

図３は、赤池情報量基準を図１と同じ模擬サンプルデータに対して適用した結果を示す図である。上記のような混合ガウシアンモデルにおける問題を回避するため、各種情報量規準を混合ガウシアンモデルに対して適用することが考えられる。赤池情報量基準は、基本的な情報量基準としてよく知られている。

赤池情報量基準を最適化において用いる場合は、多くの情報量を与えて最適化を実施するほど多くのペナルティ値を課す。これを混合ガウシアンモデルに対して適用する場合、クラスタ数を増やすとその分だけペナルティを課すので、必ずしもクラスタ数を増やすほど良いフィッティングが得られるとは限らないことになる。したがって、評価値が極値を有するクラスタ数が最適であろうと推測することができる。

ＭａｔＬａｂ２００８ａが搭載しているＥＭアルゴリズムを用いて、赤池情報量基準を混合ガウシアンモデルに対して適用した最適化計算を実施したところ、図３に示すように明確な極値を得ることはできなかった。したがって、混合ガウシアンモデル、および赤池情報量基準をこれに適用した修正手法のいずれも、最適なクラスタ数を決定することは困難であることが分かった。

その他のクラスタリングアルゴリズムとしては、サポートベクトルマシン法、ｋ−ｍｅａｎｓ法などがある。しかしこれらの手法は、事前にクラスタ数が分からないデータに対して適用しても、有効な結果を得ることは困難である。仮に、これら手法によって最適なクラスタ数が得られたとしても、そのクラスタリングの信頼度を数値によって評価することはやはり困難である。

事前に情報を与える必要がないデータ解析手法として、多くのデータマイニング手法が知られている。たとえば、非特許文献２において採用されている自己組織化マップ（Ｓｅｌｆ Ｏｒｇａｎｉｚｉｎｇ Ｍａｐｓ）がある。しかし、自己組織化マップを用いてクラスタリングを実施したとしても、クラスタリング結果の信頼度を得ることができない。

上記課題の他にも、単一細胞解析において得られるサンプルデータは実験誤差が無視できないという課題がある。実験誤差を多く含むサンプルデータは、誤差の少ないサンプルデータのクラスタリング結果から外れてしまうため、いずれのクラスタに属するのか、あるいはそもそも独立のクラスタとみなすべきであるのかを判断することが難しい。したがって、実験誤差を含むサンプルデータについて、有意な分解能でクラスタリングを実行することも、細胞解析／診断装置として重要であると考えられる。

＜実施の形態１＞
以上、従来のデータ解析手法およびその課題について説明した。以下では本発明の実施形態１に係るデータ解析装置について説明する。単一細胞ごとの生体分子に関する解析、特に遺伝子発現解析において得られるデータは、各細胞の生体分子の定量値を要素値とする行列を用いて表現することができる。個々の細胞に対する各遺伝子の発現量のサンプルデータは、遺伝子数をｎ、細胞数をｍとすると、ｍ行ｎ列の行列で表示できる。以下ではこの形式によって記述されたサンプルデータを前提とする。

図４は、本実施形態１に係るデータ解析装置がサンプルデータをクラスタリングする処理の概略フローを示す図である。以下、図４に示す各ステップを説明する。各ステップの詳細については後述の図５〜図６を用いて改めて説明する。

（図４：ステップＳ４０１〜Ｓ４０２）
データ解析装置は、図１で説明した形式のサンプルデータを取得する（Ｓ４０１）。データ解析装置は、サンプルデータを取得したときの実験誤差に関するデータを取得する（Ｓ４０２）。これらのデータは、オペレータがデータ解析装置へ入力してもよい。実験誤差に関するデータとは、サンプルデータが実験誤差に起因してどの程度ばらついているかを数値的に示すデータである。たとえば、各遺伝子データの標準偏差ベクトル値を実験誤差データとすることができる。

（図４：ステップＳ４０３）
データ解析装置は、実験誤差に起因していずれのクラスタにも属さなくなっている例外データをいずれかのクラスタに包含させるため、以後のクラスタリング処理においてクラスタ境界を緩和する。具体的には、クラスタリング空間において例外データから所定範囲離れた位置にいずれかのクラスタが存在する場合は、その例外データは当該クラスタに属するものとみなす。上記所定範囲のことを、本発明においてはＣＲ（Ｃｌｕｓｔｅｒｉｎｇ Ｒｅｓｏｌｕｔｉｏｎ）値と呼ぶことにする。例外データは実験誤差によって生じると考えられるため、ＣＲ値は実験誤差データが記述している実験誤差以上の値として設定する。例えば実験誤差データを誤差の標準偏差σによって表す場合、ＣＲ値はσ〜４σ程度の値とすることができる。

（図４：ステップＳ４０４）
データ解析装置は、ＣＲ値を掃引しながら各ＣＲ値について以下のステップＳ４０５を実施する。ＣＲ値の掃引範囲としては、実験誤差データを誤差の標準偏差σによって表す場合、例えば上述のようにσ〜４σなどとすればよい。

（図４：ステップＳ４０５）
データ解析装置は、クラスタ数ｋを２〜予想最大値まで仮設定して実際にクラスタリングを実施しながら、各クラスタリング結果についてクラスタ数の最適度を評価する。具体的には、サンプルデータが各クラスタに属する確率の対数尤度および属さない確率の対数尤度を用いて、現在のクラスタ数の尤もらしさを評価する。

（図４：ステップＳ４０６）
データ解析装置は、ステップＳ４０５において求めたクラスタ数の尤もらしさが極値を取るとき、そのクラスタ数が最適であるものとみなし、そのクラスタ数を最終的なクラスタ数として採用する。

（図４：ステップＳ４０７）
データ解析装置は、ステップＳ４０６において決定したクラスタ数に基づくクラスタリング結果とともに、そのクラスタリング結果の信頼度を出力する。クラスタリング結果の信頼度としては、ステップＳ４０５において求めたクラスタ数の尤もらしさの値を用いることができる。

図５は、ステップＳ４０５の詳細を説明する処理フロー図である。以下、図５の各ステップについて説明する。

（図５：ステップＳ５０１）
データ解析装置は、与えられた仮のクラスタ数ｋについて、サンプルデータを暫定的にクラスタリングする。本ステップにおけるクラスタリング手法は、例えば階層的クラスタリング法やｋ−ｍｅａｎｓ法など任意の手法でよい。

（図５：ステップＳ５０２）
データ解析装置は、サンプルデータをクラスタリングしたデータセットをｋ個分複製する。複製された各データセットは、それぞれ以下のステップにおいて各クラスタリング結果の尤もらしさを評価するために用いる。データ解析装置は、以下のステップにおいて用いるカウンタｉを初期化する（ｉ＝１）。

（図５：ステップＳ５０３）
データ解析装置は、ｉ番目（ｉ＝１〜ｋ）のクラスタについて、ｉ番目の複製データセットを用いて以下のステップを実施する。ｉ番目の複製データセットに関してはｉ番目のクラスタのみを残し、その他データは全てｉ番目のクラスタ外に属するものと仮定し、ｉ番目のクラスタ以外のクラスタについては削除する。すなわち、ｉ番目のクラスタに属さないデータについては、全て単一のクラスタに属するものと仮定する。

（図５：ステップＳ５０４）
データ解析装置は、ｉ番目のクラスタが例外データによって構成されているものであるか否かを判定する。例外データの例については後述の図６で説明する。例外データでなければステップＳ５０５〜Ｓ５０７を実施し、例外データであればステップＳ５０８〜Ｓ５０９を実施する。

（図５：ステップＳ５０４：例外データであるか否かの判定基準その１）
ｉ番目のクラスタがその内部に十分なサンプルデータ数を保持している場合は、そのクラスタは例外データによるものではないと判定する。この場合は、クラスタに属するサンプルデータから計算される相関行列などを用いて決定したデータ構造は、信頼度が高いと認められるからである。サンプルデータ数が十分であるか否かの閾値ｔｈは、あらかじめ定めておく。例えばクラスタ内のサンプルデータ数が２個以下であれば例外データとみなす、などの判定基準が考えられる。

（図５：ステップＳ５０４：例外データであるか否かの判定基準その２）
閾値ｔｈは、例えばある範囲内でランダムに選択することもできる。あるいは、適当な確率分布を仮定して、その確率分布を前提としてランダムに選択することもできる。この場合は、確率分布の中央に位置するクラスタ数が最も選択され易くなる。確率分布のパラメータは任意に定めることもできるし、最適化計算によって望ましい確率分布を決定することもできる。

（図５：ステップＳ５０５〜Ｓ５０７の総括）
データ解析装置は、ｉ番目のクラスタ内に属するサンプルデータの分布を用いてそのサンプルデータの確率分布を決定する。データ解析装置は、各サンプルデータがｉ番目のクラスタに属する、あるいは属さないと判定した判断が正しいか否かについて、上記確率分布を用いてその妥当性を評価する。具体的な手法は以下の通りである。

（図５：ステップＳ５０５〜Ｓ５０７）
データ解析装置は、ｉ番目のクラスタのクラスタ中心（クラスタ内に属するサンプルデータの平均）とクラスタ内のサンプルデータの標準偏差を算出し、サンプルデータを規格化する（Ｓ５０５）。データ解析装置は、ｉ番目のクラスタ内のサンプルデータの相関行列の逆行列を計算する（Ｓ５０６）。データ解析装置は、各サンプルデータとｉ番目のクラスタ中心との間のマハラノビス距離を計算する（Ｓ５０７）。クラスタ中心からのマハラノビス距離を計算する理由については後述の図６〜図７で説明する。

（図５：ステップＳ５０８〜Ｓ５０９）
データ解析装置は、ｉ番目のクラスタのクラスタ中心を算出し、クラスタ中心とＣＲ値を用いてサンプルデータを規格化する（Ｓ５０８）。データ解析装置は、各サンプルデータとｉ番目のクラスタ中心との間のユークリッド距離を計算する（Ｓ５０９）。クラスタ中心からのユークリッド距離を計算する理由については後述の図６〜図７で説明する。

（図５：ステップＳ５１０：ステップ１）
データ解析装置は、ステップＳ５０１においてｉ番目のクラスタに属さないと判定したサンプルデータについて、クラスタ中心からの距離が離れるほど当該サンプルデータがｉ番目のクラスタに属さない確率が高くなるような確率分布関数を用いて、当該サンプルデータがｉ番目のクラスタに属さない確率を算出する。同様にステップＳ５０１においてｉ番目のクラスタに属すると判定したサンプルデータについて、クラスタ中心からの距離が離れるほど当該サンプルデータがｉ番目のクラスタに属する確率が低くなるような確率分布関数を用いて、当該サンプルデータがｉ番目のクラスタに属する確率を算出する。例えば、前者については自由度ｎのｘ２乗分布の累積確率分布関数にしたがって確率値を計算し、後者については１から上記累積確率分布関数を引いた関数にしたがって確率値を計算する。この関数の例については後述の図６〜図７で例示する。

（図５：ステップＳ５１０：ステップ２）
データ解析装置は、上記関数にしたがって計算した確率値の対数尤度を計算することにより、各サンプルデータがｉ番目のクラスタに属することの尤もらしさ、または属さないことの尤もらしさを計算する。すべてのサンプルデータおよびすべてのクラスタについてこの対数尤度の和を求め、クラスタ数ｋで除算することにより、現在のクラスタ数ｋの尤もらしさを計算する。本ステップの次は、ｉ＋１番目のクラスタについてステップＳ５０３から同様の処理を実施する。

（図５：ステップＳ５１０：補足）
対数尤度の評価に用いているのはクラスタリングが妥当である確率の評価値であるため、最適パラメータが得られた時の対数尤度の値から、クラスタリングの信頼度を確率値で出力することが可能である。

（図５：ステップＳ５１１：最適化ループおよびクラスタ数掃引ループ）
データ解析装置は、例えばモンテカルロ法などのような最適化手法を用いて、ステップＳ５１０で計算した対数尤度が小さくなるようにクラスタリングを繰り返し実施することにより、最適なクラスタリング結果と最適なクラスタ数を決定する。最適化手法として例えばモンテカルロ法を用いる場合は、各クラスタに属するサンプルデータをランダムに入れ替えながらステップＳ５０３から同様の処理を実施する。最適化ループの終了条件は、例えばステップＳ５１０で計算した現在のクラスタ数ｋの尤もらしさが所定閾値に達した時点などとすればよい。最適化ループを完了した後はクラスタ数ｋをインクリメントしてステップＳ５０１に戻り、同様の処理を実施する。

（図５：ステップＳ５１１：ＣＲ値掃引ループ）
データ解析装置は、現在のＣＲ値について最適化ループおよびクラスタ数掃引ループを終了した後は、ＣＲ値をインクリメントしてステップＳ５０１に戻り、同様の処理を実施する。ＣＲ値をインクリメントする幅は、想定しているＣＲ値の最小値と最大値の差に応じて適宜定める。

図６は、ステップＳ５０５〜Ｓ５０７の処理イメージを示す図である。ここでは図１（ｂ）と同様に、２つの主因子に対応する軸を設定したクラスタリング空間を例示した。図６に示すクラスタリング例においては、クラスタ中心に比較的多くのサンプルデータが集まっているため、当該クラスタ内に含まれるデータは実験誤差の影響を受けた例外データではないと仮定する。

クラスタの形状は必ずしもクラスタリング空間において円形ではないため、図６（ａ）の左に示すデータ分布から右に示すデータ分布へと線形変換を実行し、変換後のクラスタ形状においてクラスタ中心から各サンプルデータまでの距離を求める。これは、各データのクラスタの中心からのマハラノビス距離を求めることに相当する。換言すると、クラスタがマハラノビス田口法における単位空間を形成していると仮定することに対応する。

ステップＳ５０５〜Ｓ５０７においては、ｉ番目のクラスタに属さないと仮判定したサンプルデータ（例えば図６（ａ）の左上のデータ群）については、図６（ｂ）に示すｘ^２累積確率分布を用いて、そのサンプルデータがｉ番目のクラスタに属さない確率を計算する。ｉ番目のクラスタに属すると仮判定したサンプルデータ（図６（ａ）の円内のデータ群）については、図６（ｂ）に示す１−ｘ^２累積確率分布を用いて、そのサンプルデータがｉ番目のクラスタに属する確率を計算する。これにより、ｉ番目のクラスタに属すると仮判定したサンプルデータについては、クラスタ中心に近いほど確率値が高くなり、ｉ番目のクラスタに属さないと仮判定したサンプルデータについては、クラスタ中心から離れるほど確率値が高くなる。

図７は、ステップＳ５０８〜Ｓ５０９の処理イメージを示す図である。ここでは図１（ｂ）と同様に、２つの主因子に対応する軸を設定したクラスタリング空間を例示した。図６に示すクラスタリング例においては、中央のクラスタに属するデータ数が少ないため（例外１においては２個、例外２においては１個）、例外データであると仮判定される。これはすなわち、ｉ番目のクラスタに属するデータ数が少ないため、そのクラスタ内のデータ構造を十分に決定することができない場合に相当する。

ステップＳ５０８〜Ｓ５０９においては、ＣＲ値を当該クラスタのサイズであると仮定する。ｉ番目のクラスタの周辺に存在する、同クラスタに属さないと判定されたデータ個数が所定個数未満である場合（例外１）、例外データは独立したクラスタを構成している可能性が高い。この場合、例外クラスタに属さないサンプルデータはその確率値が大きくなる。結果として、例外データが独立したクラスタであることの尤もらしさが高く評価されることになる。

一方、ｉ番目のクラスタの周辺に、同クラスタに属さないと判定されたデータ個数が所定個数以上存在する場合（例外２）、例外データは本来他のクラスタに属していたが実験誤差によってクラスタから外れた可能性が高い。この場合、近傍のクラスタが例外クラスタに属する確率値が相応に高くなる。結果として、例外データが独立したクラスタであることの尤もらしさが低く評価されることになる。

図７に示すように、本発明におけるＣＲ値は、実験誤差によってクラスタから外れたサンプルデータを元のクラスタに戻す機能を備える。そのためＣＲ値は、実験誤差よりも大きな値とし（実験誤差より小さいと誤差を修正できないため）、他の医学的、生物学的な知見から得られる制限範囲内で設定することが望ましい。

図８は、図１で説明した模擬サンプルデータについてＣＲ値を１、２、４（それぞれ４×σ、８×σ、１６×σに対応）に掃引し、クラスタ数の尤もらしさを表す対数尤度を計算した結果を示す。クラスタ数ｋ＝６において対数尤度が最小値を示し、その最小値はＣＲ値を変えても安定であることが観察される。すなわち、クラスタ数の対数尤度が極値となるクラスタ数が最適であると判定することができる。

＜実施の形態１：まとめ＞
以上のように、本実施形態１に係るデータ解析装置は、実験誤差に基づき定めた、クラスタ境界を緩和するクラスタ範囲パラメータ（ＣＲ値）を用いて、いずれのクラスタにも属さないと仮判定された例外データが他のクラスタに属するか否かを判定する。これにより、実験誤差に起因する例外データが生じる場合であっても、精度よくクラスタリングを実施することができる。

また、本実施形態１に係るデータ解析装置は、クラスタ中心からのマハラノビス距離またはユークリッド距離を基準として、ＣＲ値を掃引しながらクラスタリング結果の尤もらしさを評価し、その値が極値を取るクラスタ数を最適とみなす。これにより、最適なクラスタ数があらかじめ分かっていない場合であっても、良好なクラスタリング結果を得ることができる。

また、本実施形態１に係るデータ解析装置は、クラスタリング結果の信頼度をクラスタリング結果と併せて出力する。これにより、特定種類に属する細胞数やその種類に属する遺伝子発現マーカによる診断精度を向上させることができる。すなわち、従来は単一細胞解析以外の方法においては複数種類の細胞に対して医学的評価を実施していたところ、本実施形態１に係るデータ解析装置のクラスタリング結果に基づき特定集団に対するバイオマーカを評価することによって、バイオマーカによる疾患の種類や状態の判定の精度が向上するものと期待される。

本実施形態１に係るデータ解析装置によって決定したクラスタは、それぞれサンプルデータから導かれた細胞の種類に対応する。最適なクラスタ数を決定できるデータ解析装置は、ひとつの細胞種と別の細胞種との間の境界を明示できる能力を有すると考えられる。それゆえ、細胞の解析や診断で得られたサンプルデータから細胞種毎の個数を出力することができる。また、個々の細胞種において、特定の遺伝子発現量などに代表されるバイオマーカの平均値などの統計量を出力するための母集団を決定することもできる。さらに、その決定には信頼度が付与されているため、ある信頼度以下のデータを採用しないなどの判定も容易である。

本実施形態１に係るデータ解析装置は、具体的には、個々の細胞の特性を１細胞ずつ解析するバイオ・医療分野の解析・診断装置において適用することができる。特に、血球を解析対象とする解析／診断装置、または尿中の細胞を対象とする解析／診断装置、または組織切片を対象とする解析／診断装置において適用することができる。以下の実施形態についても同様である。

＜実施の形態２：装置構成＞
本発明の実施形態２では、実施形態１で説明したデータ解析装置の具体的な適用例として、単一細胞中の遺伝子発現解析によるデータによって細胞を分類するデータ解析装置について説明する。

本実施形態２においては、単一細胞中の遺伝子発現解析を実施するため、非特許文献１に記載されている方法に基づき、磁性ビーズ上に単一細胞中のｃＤＮＡライブラリを作製し、単一細胞中の０．５ｐｇという微量なｍＲＮＡをｑＰＣＲ装置によって定量する方法を用いた。

図９は、本実施形態２に係るデータ解析装置９０６の機能ブロック図である。データ解析装置９０６は、演算部９０４、データ入出力部９０５、サンプルデータ入力部９０８、実験誤差データ入力部９０９を備える。各データを測定する測定系についても図９内に併記した。

図９は、実験誤差データを取得する動作と、クラスタリング対象となる単一細胞解析とを平行して実施する場合における測定系の構成を示した。機能ブロック９０１は、クラスタリング対象となるサンプルデータを取得する機能部である。機能ブロック９０２は、実験誤差を取得する機能部である。機能ブロック９０３は、両者に共通する処理を実施する機能部である。正確な実験誤差を評価するためには、実験誤差を評価するためのデータとクラスタリング対象とするサンプルデータは共通である部分が多い方が望ましい。実験誤差に関する全部または一部のデータを事前に取得している場合は、そのデータをデータベース９０７に保存しておき、クラスタリングを実施するとき演算部９０４がこれを読み出すようにしてもよい。

機能ブロック９０１は、細胞を１つずつ採取し、個々の細胞ごとに反応容器に分注し、この細胞を分注した反応容器に、ｍＲＮＡを抽出するためのポリＴプローブ付きビーズとリシスバッファを導入する。

機能ブロック９０２は、多数の細胞を反応容器に導入した後、機能ブロック９０１と同様にｍＲＮＡ抽出を実施し、その後ｍＲＮＡ溶液を希釈し、単一細胞相当量を採取して、ポリＴプローブ付ビーズ溶液にこの溶液を分注する。

機能ブロック９０３は、リシスバッファを除去し、逆転写酵素を含む反応溶液を分注し、逆転写反応後に反応溶液を除去し、ｍＲＮＡの分解酵素を加えて洗浄後、ｑＰＣＲ試薬を導入し、ＰＣＲのための温度サイクルを印加しながら蛍光測定を実施することによって定量を実施する。

個々のサンプル中の遺伝子発現量は、蛍光強度がある閾値を超えたサイクル数Ｃｔ値に基づき算出される。分子数が既知の複数のＤＮＡに対してｑＰＣＲ定量を実行することによって、Ｃｔ値は分子数に読み替えられる。実験誤差は、定量を実施するまでにおける単一細胞のサンプリング以外のすべてのプロセス内の誤差を含む。定量値はガウス分布に従っているほうが望ましいので、分子数の対数をとった値をサンプルデータとして演算部９０４に出力する。実験誤差データについても同様である。

サンプルデータ入力部９０８は、機能ブロック９０１、９０３を介してサンプルデータを取得する。実験誤差データ入力部９０９は、機能ブロック９０２、９０３を介して実験誤差データを取得する。演算部９０４は、これらのデータを用いて、実施形態１で説明した処理フローにしたがってクラスタリングを実施する。データ入出力部９０５は、これらデータの入出力を制御する。

本実施形態２においては、仮クラスタ内の要素数が３より小さい（２以下）場合に、ステップＳ５０８〜Ｓ５０９を実施することとした。生物学的あるいは医学的に特定の遺伝子の発現量（サンプルデータの値）が変動すると知られている場合、その変動幅をＣＲ値として用いてもよい。このような、遺伝子とＣＲ値の対応関係やサンプルデータのカテゴリとＣＲ値の対応関係をデータベース９０７に保存しておき、データ解析装置９０６に入力された遺伝子やサンプルデータのカテゴリのなかでデータベース９０７に保存しておいたものと合致する該当するものがあれば、演算部９０６はそれをデータベース９０７から読み出す。あるいは、実験誤差データに代えてステップＳ５０４における判定閾値をデータベース９０７にあらかじめ格納しておき、データ解析装置９０４に入力されたサンプルデータや遺伝子情報に基づいて適当な値を用いるようにしてもよい。

本実施形態２では、ステップＳ５０４における判定閾値は、クラスタ内のサンプルデータ数に対応する固定した値としたが、判定閾値を複数設けておき、各閾値をそれぞれ用いてクラスタリングを実施し、もっとも対数尤度の低い（最も尤もらしい）閾値を選択するようにしてもよい。さらには、閾値をランダムに選択するか、あるいは閾値がある確率分布にしたがう前提の下でランダムに選択してもよい。このときの確率分布関数はあらかじめ複数種類をデータベース９０７に格納しておき、サンプルデータの内容に応じて適当な確率分布関数を選択すればよい。

＜実施の形態２：サンプリング結果＞
以下では、本実施形態２に係るデータ解析装置９０６によるサンプリング結果の例を説明する。サンプルとしてマウスの間葉系幹細胞（Ｃ３Ｈ１０Ｔ１／２）を９２個用いた。すなわち、細胞を取得したマウスの骨の分化誘導状態を知るためにこのような細胞を採取した。もちろん、他の目的のための他の生体（人を含む）から、別の細胞（たとえば、血中の免疫細胞やがんの組織切片など）を採取してもよい。サンプルデータとして用いる遺伝子発現データは、Ｂｇｌａｐ１、Ｃｏｌ１ａ１、Ｐｐａｒｇ、Ｃｏｌ２ａ１、Ｅｅｆ１ｇの５種類の遺伝子についての定量値とした。遺伝子の種類およびその数は１例であり、可能性のある遺伝子すべてを計測対象にしてもかまわない。

次に、実験誤差を評価するため、５×１０^５個程度以上の多数の細胞から抽出したｍＲＮＡを単一細胞相当量である０．５ｐｇずつ９６個分サンプリングし、５×１０^５個の細胞の遺伝子発現量の平均値を計測するためｑＰＣＲ装置を用いて定量した。定量値のばらつきはｃＤＮＡライブラリを作製するときのハンドリング誤差や、ｑＰＣＲ定量時の定量誤差から構成される実験誤差の総量である。この誤差は遺伝子ごとに異なる。この誤差を５次元ベクトルとして評価した。

実験誤差は以下のような方法によって取得することが望ましい。まず、多数の細胞から抽出したｍＲＮＡサンプルを作成するとき、細胞数の異なる複数のサンプルを準備し、これらのサンプルから単一細胞相当量のｍＲＮＡを採取し、ｃＤＮＡを構築後、各細胞数および各遺伝子に関するｑＰＣＲ定量時の誤差を定量する。そして、細胞数無限大の場合の値を外挿によって推定する。実際には、細胞数の逆数を横軸に実験誤差（標準偏差）を縦軸にとって、ｙ切片の値を実験誤差の推定値とする。

得られた実験誤差に基づき、許容できるＣＲ値を定める。このＣＲ値（ベクトル値）をσとする。本実施形態２においては、実験誤差と最小ＣＲ値が一致するものとする。ただし、特定の遺伝子の発現は、同じ細胞状態であっても時間的に変化しており、実験誤差によるデータのばらつき以外に生物学的にもデータの変動幅が分かっている場合がある。これらの変動をクラスタリングにおいて用いたくない場合は、σの値をサンプルデータに基づいて一部変更しても構わない。具体的には、たとえば、特定の遺伝子について１０倍程度の変動がある場合、この遺伝子のみについて、実験誤差ではなくこの変動幅の数値をＣＲ値として用いてもよい。ただし、実験誤差が生物学的あるいは医学的変動より十分小さい場合に限られる。

図１０は、本実施形態２において、クラスタ数の尤もらしさを表す対数尤度を計算した結果を示す。実験誤差をσとして、ＣＲ値は１σから２．１σまでについて評価した。すべてのＣＲ値について、クラスタ数ｋ＝１５のとき、安定的に極値が得られた。

図１１は、階層的クラスタリングの結果と本実施形態２に係るデータ解析装置９０６によるクラスタリング結果を比較する図である。１５個のクラスタは階層的クラスタリング結果に対応するトーナメント図の下に四角の箱で表現している。１５個のクラスタは、５個の有意なクラスタと１０個の例外クラスタから構成されることが分かった。５個の有意なクラスタは階層的クラスタリングの結果と一致する。本実施形態２により、有意なクラスタと例外クラスタの境界が明確になったことが分かる。

図１２は、各クラスタ（１〜５の番号で表示）のＢｇｌａｐ１、ＰＰａｒｇ, Ｃｏｌ２ａ１の遺伝子発現量を示す図である。縦軸は、単一細胞中の各遺伝子の分子数である。横軸はクラスタ番号を示している。黒いバーは分化誘導剤あり、灰色のバーは分化誘導剤なしに対応する細胞についての遺伝子発現量である。分化誘導剤あり・なしの区別は事前に分かっている細胞情報の例である。

図１２に示すように、分化誘導剤あり・なしは、個々のクラスタのなかでは細胞を区別するための重要な情報でないことが分かる。一方、３つの遺伝子の発現プロファイル、すなわち、（Ｂｇｌａｐ１，Ｐｐａｒｇ，Ｃｏｌ２ａ１）の発現量を＋−で表現すると、クラスタ１は（−，＋，＋）、クラスタ２は（＋，＋，＋）、クラスタ３は（−，−，＋）、クラスタ４は（−，＋，−）、クラスタ５は（＋，−，−）という明確な特徴を持っていることがわかる。同時にクラスタごとの遺伝子発現量の変動量は、上記特徴を識別するには十分でないことがわかる。事前情報との対応は各クラスタに属する細胞数に反映されるので、その意味においても正確なクラスタリングが重要であることがわかる。

図１３は、図１２で説明した３つの遺伝子について、クラスタリングを実施せず、分化誘導ありの場合となしの場合について発現量を示した図である。図１３に示すようにＢｇｌａｐ１については事前情報に対応した変化が見られるが、ＰｐａｒｇおよびＣｏｌ２ａ１については事前情報に対応した発現量の変化がなく、これら２つの遺伝子に関しては情報が得られない。以上から、クラスタリングによって、より詳細な生体の遺伝子発現に関する情報が得られることが分かる。

以上のように、本実施形態２によれば、個々の細胞中の遺伝子発現量を計測し、クラスタリングによって、細胞を構成する生体の情報を得ることができる。すなわち本実施形態２に係るデータ解析装置９０６は、生体中にどのような種類（クラスタ）に属する細胞がどの程度の数だけ存在するかを推定することによって生体の状態を推定するための装置である。計測対象となっている生体の健康状態が変化すると、クラスタの種類とクラスタに属する細胞の数が変化する場合において、本実施形態２は有効である。

＜実施の形態３＞
図１４は、本発明の実施形態３に係るデータ解析装置９０６の機能ブロック図である。本実施形態３では、遺伝子発現の定量方法として、大規模ＤＮＡシーケンサを用いる構成例を説明する。これにより、計測対象とする遺伝子数を大幅に増やすことができる。

機能ブロック９０１は、クラスタリング対象サンプルをプレート上のポリＴプローブ付ビーズの入った反応容器に１つずつ分注し、反応容器内で細胞を破砕し、ｍＲＮＡをビーズ表面でトラップすることによって抽出する。

機能ブロック９０２は、実験誤差測定のため、複数の遺伝子について既知量のｍＲＮＡを個々の反応容器に分注する。既知量ＲＮＡを導入した、反応容器に対応するシーケンシングデータからその標準変偏差を算出することによって、サンプル処理とシーケンシングに関する実験誤差を定量することができる。

機能ブロック９０３は、逆転写反応とｍＲＮＡ分解を実施する。このとき、ｃＤＮＡの末端に一括増幅用プライマを導入し、次にこのプライマを用いて一括増幅を実施する。さらに、断片化し、各反応槽ごとに配列のことなる細胞認識タグをもった増幅プライマ（プライマ配列はすべての容器で共通）を用いて一括増幅を実施し、シーケンシングライブラリを構築する。シーケンシングライブラリの末端には容器ごとすなわち細胞ごとに異なる配列のタグが挿入されているので、以下の処理においてはサンプルを混合することができる。混合したサンプルを大規模シーケンサでシーケンシングするために、エマルジョンＰＣＲやブリッジＰＣＲなどの個別増幅を用いてシーケンシングを実行する。シーケンシングは、蛍光計測を用いる装置、ＦＥＴを用いる装置、ナノポアを用いる装置などいずれの装置を用いてもよい。断片化したサンプルから得られたシーケンシングデータを既知の遺伝子配列にマッピングすることによって、どの遺伝子のどの領域の配列が計測できたかを判定する。その後、データを遺伝子毎に集計し、遺伝子毎の発現量データを算出する。このときの算出アルゴリズムは当該専門家が一般的に用いるアルゴリズムを用いてよい。その結果、細胞毎・遺伝子毎の発現量データが得られる。

クラスタリング手順については実施形態１〜２と同様であるが、計測遺伝子数が数万程度と非常に大きいため、細胞をより詳細に分類することができる。

同様の処理を用いて、各細胞のゲノムデータを解析することもできる。データ解析装置９０６に入力されるデータとしては、ゲノムを複数の領域に分割し、各分割領域についてカウントした変異数である。この場合の計測目的は、例えば癌の組織切片からの細胞を計測することによって、癌の発生や転移のメカニズムを解明すること、あるいは分子標的薬を選択するための診断、などが考えられる。

ゲノムデータを解析することによって得られるデータについて補足する。ゲノム全体またはその一部を単一細胞ごとに配列解析し、（ｍＲＮＡの配列解析におけるデータを用いることを想定）たとえば５０ｋ塩基ごとに領域を分けて、変異を計測する。計測対象は、たとえば一塩基置換、欠失、挿入、遺伝子コピー数異常、などである。入力データはそれぞれの変異数である。実験誤差は、変異がないサンプルを人工的に作製し、それを評価することによって、評価することができる。

ゲノムを直接配列決定するためには、ｍＲＮＡ抽出の替わりにＲＮＡ分解、ＤＮＡ抽出後、断片化、ポリＡテイリングを酵素試薬の添加によって実行する必要がある。その後はｍＲＮＡの処理と同様である。

＜実施の形態４＞
本発明の実施形態４では、細胞中の遺伝子の発現量（分子数）を定量することによって細胞を特徴付ける代わりに、免疫染色法などで得られた細胞サンプルを蛍光顕微鏡でイメージングして、蛍光強度と分子数の対応データまたは単分子蛍光カウンティングによって細胞中の蛋白質量を定量することにより細胞を分類する構成例を説明する。

本実施形態４においては、遺伝子の種類が蛋白の種類に対応する。細胞ごとの蛋白量（分子数）を、サンプルデータとしてデータ解析装置９０６に入力する。免疫染色などのサンプル作製時は、蛍光計測時の誤差を評価して実験誤差としてデータ解析装置９０６に入力する。これにより、細胞のクラスタリングを実行することができる。

図１５は、本実施形態４に係るデータ解析装置９０６の機能ブロック図である。機能ブロック９０１は、細胞サンプルとして組織切片を採取し、脱パラフィン処理後に免疫染色を実施する。このとき、免疫染色と異なる蛍光波長の既知量の蛍光色素を細胞内に導入すると同時に核染色と細胞膜染色を実施する。７色の多色の蛍光計測を実施することによって、ターゲット蛋白（３種類）と実験誤差計測用色素、核、細胞膜に対応する像を得る。ターゲット蛋白の種類を増やす場合には計測色を増やす必要がある。細胞膜と核を認識し、細胞膜中に核が１つある物体を細胞と認識する。個々の細胞と認識された領域のターゲット蛋白に相当する蛍光色の強度、または単分子蛍光時の分子数を定量値として記録する。同時に実験誤差用色素の定量も実行する。実験誤差は、認識した細胞の面積あたりの強度の標準偏差を求め、細胞の平均面積から算出する。

本発明は上記した実施形態に限定されるものではなく、様々な変形例が含まれる。上記実施形態は本発明を分かりやすく説明するために詳細に説明したものであり、必ずしも説明した全ての構成を備えるものに限定されるものではない。また、ある実施形態の構成の一部を他の実施形態の構成に置き換えることもできる。また、ある実施形態の構成に他の実施形態の構成を加えることもできる。また、各実施形態の構成の一部について、他の構成を追加・削除・置換することもできる。

上記各構成、機能、処理部、処理手段等は、それらの一部や全部を、例えば集積回路で設計する等によりハードウェアで実現してもよい。また、上記の各構成、機能等は、プロセッサがそれぞれの機能を実現するプログラムを解釈し、実行することによりソフトウェアで実現してもよい。各機能を実現するプログラム、テーブル、ファイル等の情報は、メモリ、ハードディスク、ＳＳＤ（Ｓｏｌｉｄ Ｓｔａｔｅ Ｄｒｉｖｅ）等の記録装置、ＩＣカード、ＳＤカード、ＤＶＤ等の記録媒体に格納することができる。

９０６：データ解析装置、９０４：演算部、９０５：データ入出力部、９０８：サンプルデータ入力部、９０９：実験誤差データ入力部。

Claims (10)

1. 複数のサンプルデータをクラスタリング解析する装置であって、
   複数のサンプルデータを受け取るサンプルデータ入力部と、
   前記サンプルデータの実験誤差についての情報を記述した実験誤差データを受け取る実験誤差データ入力部と、
   前記複数のサンプルデータをクラスタリング空間内においてクラスタリングする演算部と、
   前記クラスタリングの結果を出力する出力部と、
   を備え、
   前記演算部は、
   前記実験誤差データが記述している前記実験誤差の範囲に応じて、前記クラスタリングを実施する際におけるクラスタ境界を緩和するためのクラスタ範囲パラメータをあらかじめ取得しておき、
   仮設定したクラスタ総数に準じて前記複数のサンプルデータをクラスタリングし、
   前記複数のサンプルデータのうちいずれのクラスタにも属さない例外データについては、
   前記クラスタリング空間上で前記例外データからさらに前記クラスタ範囲パラメータによって定められる距離だけ離れた領域がいずれかのクラスタに含まれる場合は、前記例外データがそのクラスタに属するものと判定し、前記領域がいずれのクラスタにも含まれない場合は、前記例外データが独立したクラスタを構成するものと判定し、
   前記演算部は、
   各前記サンプルデータが前記クラスタリングの結果得られた各クラスタに属する確率の尤もらしさを表す第１対数尤度と、各前記サンプルデータが前記クラスタリングの結果得られた各クラスタに属さない確率の尤もらしさを表す第２対数尤度とをそれぞれ算出する処理を、前記第１対数尤度と前記第２対数尤度を用いて算出される前記クラスタリング結果の尤もらしさが所定閾値に達するまで繰り返すことにより、最適なクラスタ総数を求め、
   得られた最適なクラスタ総数に準じて、前記複数のサンプルデータの最終的なクラスタリング結果を決定する
   ことを特徴とするデータ解析装置。
2. 前記演算部は、
   前記仮設定したクラスタ総数に準じて前記複数のサンプルデータをクラスタリングする過程において、前記サンプルデータが仮設定したクラスタに属すると仮定して前記第１対数尤度と前記第２対数尤度を算出し、
   前記仮設定したクラスタに属すると仮定した前記サンプルデータについては、前記クラスタリング空間上において、前記仮設定したクラスタの中心からの距離が離れるにしたがって、前記仮設定したクラスタに属する確率を低く評価し、
   前記仮設定したクラスタに属さないと仮定した前記サンプルデータについては、前記クラスタリング空間上において、前記仮設定したクラスタの中心からの距離が離れるにしたがって、前記仮設定したクラスタに属さない確率を高く評価する
   ことを特徴とする請求項１記載のデータ解析装置。
3. 前記演算部は、
   クラスタ内に属する前記サンプルデータの個数が所定個数以上であるか否かを基準として、前記サンプルデータが前記例外データであるか否かを判定する
   ことを特徴とする請求項１記載のデータ解析装置。
4. 前記演算部は、前記所定個数をランダムに選択する
   ことを特徴とする請求項３記載のデータ解析装置。
5. 前記演算部は、所定の確率分布にしたがって前記所定個数をランダムに選択する
   ことを特徴とする請求項３記載のデータ解析装置。
6. 前記演算部は、
   前記クラスタ範囲パラメータを掃引し、各前記クラスタ範囲パラメータの値を採用した場合における前記クラスタリングによって得られたクラスタ総数を取得し、
   前記第１対数尤度と前記第２対数尤度に基づき算出した前記クラスタリング結果の尤もらしさの値が極値となるクラスタ総数を、最適なクラスタ総数として採用する
   ことを特徴とする請求項１記載のデータ解析装置。
7. 前記演算部は、
   前記クラスタリングを実施する過程において得られる情報を用いて、前記クラスタリングの結果の信頼度指標を算出し、
   前記出力部は、
   前記信頼度指標を前記クラスタリングの結果と併せて出力する
   ことを特徴とする請求項１記載のデータ解析装置。
8. 前記演算部は、
   前記第１対数尤度と前記第２対数尤度に基づき算出した前記クラスタリング結果の尤もらしさの値を、前記クラスタリングの結果の信頼度指標として算出する
   ことを特徴とする請求項７記載のデータ解析装置。
9. 前記サンプルデータ入力部と前記実験誤差データ入力部は、細胞の解析結果に関するデータをそれぞれ前記サンプルデータおよび前記実験誤差データとして受け取り、
   前記演算部は、前記クラスタリングによって前記細胞をグループ化する
   ことを特徴とする請求項１記載のデータ解析装置。
10. 複数のサンプルデータをクラスタリング解析する方法をコンピュータに実行させるデータ解析プログラムであって、前記コンピュータに、
    複数のサンプルデータを受け取るサンプルデータ入力ステップ、
    前記サンプルデータの実験誤差についての情報を記述した実験誤差データを受け取る実験誤差データ入力ステップ、
    前記複数のサンプルデータをクラスタリング空間内においてクラスタリングする演算ステップ、
    前記クラスタリングの結果を出力する出力ステップ、
    を実行させ、
    前記演算ステップにおいては、前記コンピュータに、
    前記実験誤差データが記述している前記実験誤差の範囲に応じて、前記クラスタリングを実施する際におけるクラスタ境界を緩和するためのクラスタ範囲パラメータをあらかじめ取得するステップ、
    仮設定したクラスタ総数に準じて前記複数のサンプルデータをクラスタリングするステップ、
    前記複数のサンプルデータのうちいずれのクラスタにも属さない例外データについては、
    前記クラスタリング空間上で前記例外データからさらに前記クラスタ範囲パラメータによって定められる距離だけ離れた領域がいずれかのクラスタに含まれる場合は、前記例外データがそのクラスタに属するものと判定し、前記領域がいずれのクラスタにも含まれない場合は、前記例外データが独立したクラスタを構成するものと判定するステップ、
    を実行させ、
    前記演算ステップにおいては、前記コンピュータにさらに、
    各前記サンプルデータが前記クラスタリングの結果得られた各クラスタに属する確率の尤もらしさを表す第１対数尤度と、各前記サンプルデータが前記クラスタリングの結果得られた各クラスタに属さない確率の尤もらしさを表す第２対数尤度とをそれぞれ算出する処理を、前記第１対数尤度と前記第２対数尤度を用いて算出される前記クラスタリング結果の尤もらしさが所定閾値に達するまで繰り返すことにより、最適なクラスタ総数を求めるステップ、
    得られた最適なクラスタ総数に準じて、前記複数のサンプルデータの最終的なクラスタリング結果を決定するステップ、
    を実行させる
    ことを特徴とするデータ解析プログラム。

Images