特許 7427583 固有分子識別子のＰＣＲ後度数からのＰＣＲ前フラグメント数の推定

(19)【発行国】日本国特許庁(JP)

(12)【公報種別】特許公報(B2)

(11)【特許番号】

(24)【登録日】2024-01-26

(45)【発行日】2024-02-05

(54)【発明の名称】固有分子識別子のＰＣＲ後度数からのＰＣＲ前フラグメント数の推定

(51)【国際特許分類】

   G16B 20/10 20190101AFI20240129BHJP

   C12Q 1/686 20180101ALN20240129BHJP

【ＦＩ】

G16B20/10

C12Q1/686 Z

【請求項の数】 22

(21)【出願番号】P 2020516693

(86)(22)【出願日】2018-09-21

(65)【公表番号】

(43)【公表日】2021-01-21

(86)【国際出願番号】 EP2018075606

(87)【国際公開番号】W WO2019057895

(87)【国際公開日】2019-03-28

【審査請求日】2021-09-09

(31)【優先権主張番号】17192640.5

(32)【優先日】2017-09-22

(33)【優先権主張国・地域又は機関】EP

(73)【特許権者】

【識別番号】508160196

【氏名又は名称】レクソジェン・ゲゼルシャフト・ミット・ベシュレンクテル・ハフツング

【氏名又は名称原語表記】ＬＥＸＯＧＥＮ ＧｍｂＨ

(74)【代理人】

【識別番号】100099759

【弁理士】

【氏名又は名称】青木 篤

(74)【代理人】

【識別番号】100123582

【弁理士】

【氏名又は名称】三橋 真二

(74)【代理人】

【識別番号】100117019

【弁理士】

【氏名又は名称】渡辺 陽一

(74)【代理人】

【識別番号】100141977

【弁理士】

【氏名又は名称】中島 勝

(74)【代理人】

【識別番号】100150810

【弁理士】

【氏名又は名称】武居 良太郎

(74)【代理人】

【識別番号】100182730

【弁理士】

【氏名又は名称】大島 浩明

(72)【発明者】

【氏名】アンドレアス テュルク

(72)【発明者】

【氏名】ミヒャエル モルダシュル

【審査官】渡邉 加寿磨

(56)【参考文献】

【文献】国際公開第２０１７／０５１３８７（ＷＯ，Ａ１）

【文献】特表２００６－５０５８３２（ＪＰ，Ａ）

(58)【調査した分野】(Int.Cl.，ＤＢ名)

Ｇ１６Ｂ ５／００－９９／００

Ｇ０６Ｑ １０／００－９９／００

Ｇ１６Ｚ ９９／００

Ｃ１２Ｑ １／６８６

(57)【特許請求の範囲】

【請求項1】

サンプル中の核酸コピー数を推定する方法であって、
a) 決定すべきコピー数の核酸を含むサンプルを提供し；
b) 前記核酸に可変的ヌクレオチド標識である標識を取り付け；
c) 前記標識付き核酸を、核酸複製手順を用いて増幅させ；
d) 増幅された標識付き核酸コピーの量を決定し、各々の量が異なる標識付き核酸コピーについて決定され；
e) ステップd)の決定量に基づいて、ステップa)のサンプル中の核酸コピー数の推定量を提供するため、コンピューターシステムによりプログラムを実行し、当該システムは、
サンプル中の少なくとも２の異なる核酸種について、前記標識が１つの核酸に結合する確率を概算するプログラムを実行し；これは増幅後の異なる核酸種について前記標識の量を平均し；そして反復的にまたは段階的に
（Ａ）ステップd)において検出された異なる標識の数および１つの核酸コピーに結合する標識の確率に従った１つの核酸に結合する異なる標識の数の期待値、およびサンプル中の核酸コピーの期待数または推定増幅効率の事前反復値または既定値に基づいて、サンプル中の核酸コピーの期待数を精密化し；
または
（Ｂ）(i) ステップd)の決定量、１つの核酸に結合する標識の確率、推定増幅効率、またはサンプル中の核酸コピーの期待数に基づいて、そして前記推定増幅効率に依存した複製サイクルあたりの核酸増幅手順での増幅された標識付き核酸コピーの推定複製率に基づいて、増幅された所定の標識付き核酸コピーの量の確率分布をモデル化し、
(ii) 標識付き核酸コピーの決定量が前記標識付き核酸コピーのモデル化された確率分布に従って起こる尤度を決定し、
(iii) 推定増幅効率またはサンプル中の核酸コピーの期待数を変化させることにより、ステップ(ii)の尤度を最大化し、
(iv) 前記最大化された尤度に従ってまたは前記最大化された尤度における推定増幅効率に従って、サンプル中の核酸コピー数の推定値を提供する
各ステップを含む方法。

【請求項2】

ステップe)において実行される前記プログラムが、少なくとも４の異なる核酸種について、前記標識が１つの核酸に結合する確率を概算する、請求項１に記載の方法。

【請求項3】

ステップe)において実行される前記プログラムが、少なくとも１０の異なる核酸種について、前記標識が１つの核酸に結合する確率を概算する、請求項１に記載の方法。

【請求項4】

ステップＡ）において、サンプル中の核酸コピーの期待数の初回反復の既定値に基づいて、推定増幅効率がステップd)の異なる検出標識の数の数値から選択されるかまたは１から増幅核酸コピーの決定量までの範囲の整数より選択される、請求項１に記載の方法。

【請求項5】

ステップＢ）において、標識付きの増幅核酸コピーの推定複製率が、ガウス分布、負の二項分布、ガンマ分布、ディラックのデルタ分布もしくはガルトン・ワトソン分布の複製確率関数またはその混合形態でモデル化される、請求項１に記載の方法。

【請求項6】

１つの核酸に結合する標識の確率が
式(12)

【数1】

式(45)

【数2】

または式(46)

【数3】

に従って決定され；そして／または
推定増幅効率またはサンプル中の核酸コピーの期待数が式(27)

【数4】

に従った変換により互換的に用いられ；そして／または
ステップＡ）において、増幅された所定の標識付き核酸コピーの量の確率分布が、式(5)

【数5】

に従って決定され；または
ステップＢ）において、増幅された所定の標識付き核酸コピーの量の確率分布が、式(25)

【数6】

に従って決定され、前記数式において、
・Kは核酸コピー数を意味し；
・ｂは標識を意味し；
・Pbは１つの核酸コピーに結合する標識の確率を意味し；
・Pdは推定増幅効率を意味し；
・Ｎは増幅された核酸コピーの決定された量を意味し；
・_estは推定値であることを意味し；
・Ki又はKi+1はi又はKi+1回目の反復を意味し；
・DUはステップd)において検出された別個の標識の数を意味し；
・cは複製サイクル数を意味し；

【数7】

はPCR前の完全UMI分布を意味し、UMIは標識の一例である固有分子識別子を意味し；
・Fはプールの総数を意味し

【数8】

はプール標識を意味し、ここでプールは標識の無い同じ種の全ての核酸分子を意味し；

【数9】

は全UMIにおけるPCR後カウントの完全セットを意味し；
・kは増幅前の標識を有する核酸コピーの数を意味し；
・nは標識を有する増幅された核酸コピーの量を意味し；
・bは標識番号を意味し；
・BはUMI即ち標識の総数を意味し；

【数10】

は確率の平均を意味し；

【数11】

は

【数12】

において観測された別個のUMI即ち標識の数を意味し；

【数13】

はKのi回目の反復における式(2)

【数14】

により定義される；
請求項１～５のいずれか一項に記載の方法。

【請求項7】

ステップＡ）において、増幅された所定の標識付き核酸コピーの量の確率分布が、式(4)

【数15】

と式(5)

【数16】

に従って決定される、請求項６に記載の方法。

【請求項8】

前記核酸増幅手順において前記標識付き核酸を増幅させる効率が、100％ではない、請求項１～７のいずれか一項に記載の方法。

【請求項9】

前記核酸増幅手順において前記標識付き核酸を増幅させる効率が、核酸コピーの異なる種ごとに異なる、請求項８に記載の方法。

【請求項10】

少なくとも20の、異なる標識が前記核酸に結合される、請求項１～８のいずれか一項に記載の方法。

【請求項11】

少なくとも30の、異なる標識が前記核酸に結合される、請求項１０に記載の方法。

【請求項12】

少なくとも40の、異なる標識が前記核酸に結合される、請求項１０に記載の方法。

【請求項13】

前記核酸がDNAまたはRNAである、請求項１～１２のいずれか一項に記載の方法。

【請求項14】

前記増幅がPCRによって行われる、請求項１～１３のいずれか一項に記載の方法。

【請求項15】

前記増幅が少なくとも１の複製サイクル間行われる、請求項１～１４のいずれか一項に記載の方法。

【請求項16】

前記増幅が２，３，４，５，６，７，８またはそれ以上の複製サイクル間行われる、請求項１５に記載の方法。

【請求項17】

前記モデル化が確率分岐モデルを使用することを含む、請求項１～１６のいずれか一項に記載の方法。

【請求項18】

更なる反復が少なくとも２回前の反復に基づき、そしてその更なる反復のサンプル中の核酸コピーの期待数または推定増幅効率が区間によって選択され、ここで前記区間は、前記標識付き核酸コピーの決定量が前記区間内の標識付き核酸コピーのモデル化確率分布に従って出現する尤度を含み、そして前記標識付き核酸コピーの決定量が前記標識付き核酸コピーのモデル化確率分布に従って出現する最大尤度が、前記区間の境界線である、請求項１～１７のいずれか一項に記載の方法。

【請求項19】

サンプル中の核酸コピーの量の推定値を計算するためにコンピューターシステムによって実行されるコンピュータープログラムであって、当該計算が、請求項１のステップe)に従って標識付き核酸の増幅後の決定量に基づいて実行される、コンピュータープログラム。

【請求項20】

ステップe)の結果を可読媒体上に提示または表示させるのに適合した、請求項１９に記載のコンピュータープログラム。

【請求項21】

請求項１に定義されるステップe)をコンピューター上で実行する方法。

【請求項22】

ステップe)の結果を可読媒体上に提示または表示させることを含む、請求項１～１８および２１のいずれか一項に記載の方法。

【発明の詳細な説明】

【技術分野】

【0001】

本発明は、バイアス（偏り）のある核酸量測定値を補正する方法の分野であって、前記バイアスが増幅法により導入される分野に関する。

【背景技術】

【0002】

ポリメラーゼ連鎖反応（PCR）は変性、アニーリングおよび重合から成る多重サイクルを適用することにより核酸を増幅する。各ステップの効率は入力配列に依存するため、PCR増幅は配列特異的バイアス（偏り）を導入する可能性があり、これは異なる配列を有する異なる核酸が異なる割合で増幅されてしまうことを意味する。例えば、RNA-Seq（RNAシークエンシング、特にフラグメントまたはショットガンシークエンシング）の場合、これは不正確な遺伝子およびイソ型の定量化を引き起こし得る（図１）。

【0003】

国際公開WO 2017/051387 A1パンフレットは、遺伝情報を使っておよび遺伝情報エントロピーを観測することによって核酸配列の複雑性を決定する方法に関する。

【0004】

固有分子識別子（UMIs）のような、固有分子バーコードとも称される標識は、PCR複製物を同定しかつ配列特異的なPCRバイアスの影響を減らすために開発された。UMIは、PCR増幅前にDNAフラグメントに連結される各位置が、大部分ランダムヌクレオチド分布を有するオリゴヌクレオチドである。UMIが均等分布しそしてそれらの数が同一断片のプールよりも相当大きい場合、同じUMIが該プール中の２つの異なるフラグメントに連結される確率は低くなる。この場合、PCR後のプール中の異なるUMIの数は、PCR前のそれのフラグメントの数と同じである。よって、異なるプール間のPCR効率の変動は、PCR後の各プール中のフラグメント数よりもむしろ異なるUMIの数をカウントすることによって補正することができる（図２）。しかしながら、この単純なカウントルーティンは、例えば、UMIオリゴヌクレオチド内の所定位置でのUMI生成または分布におけるライゲーションの偏りや一様でないヌクレオチド挿入のため、UMIが均等分布（一様分布）しない場合に、不正確なものになる。不均等なUMI分布は、UMIセットの不十分なサイズと同じ影響を与える。どちらの場合も、同一UMIが２つの断片に対して高確率で連結され、PCR前のフラグメントコピー数の過小評価を引き起こすだろう。

【0005】

従って、増幅後の量の測定値に基づいて核酸量を推定できかつそのような増幅バイアスを補正することができる、改良された方法への必要性が存在する。

【発明の概要】

【課題を解決するための手段】

【0006】

本発明は、サンプル中の核酸コピー数を推定する方法であって（増幅後の概算を利用する）、
a) 決定すべきコピー数の核酸を含むサンプルを提供し；
b) 前記核酸に可変的標識を取り付け；
c) 標識付き前記核酸を、核酸複製手順を用いて増幅させ；
d) 増幅された標識付き核酸コピーの量を決定し、各々の量が異なる標識付き核酸コピーについて決定され；
e) ステップd)の決定量に基づいて、ステップa)のサンプル中の核酸コピー数の推定値を提供し、これは
サンプル中の少なくとも２、好ましくは少なくとも４、少なくとも１０以上の異なる核酸種について、前記標識が１つの核酸に結合する確率を概算し；この確率は、増幅後の異なる核酸種について前記標識の量を平均し；そして反復的にまたは段階的に
（Ａ）ステップd)において検出された異なる標識の数および１つの核酸に結合する異なる標識の数の期待値（１つの核酸コピーに結合する標識の確率に従う）およびサンプル中の核酸コピーの期待数または推定増幅効率の事前反復値もしくは既定値に基づいて、サンプル中の推定核酸コピー数を精密化し；
または
（Ｂ）(i) ステップd)の決定量、１つの核酸に結合する標識の確率、推定増幅効率、またはサンプル中の核酸コピーの期待数に基づいて、そして前記推定増幅効率に依存した複製サイクルあたりの核酸増幅手順において増幅された標識付き核酸コピーの推定複製率に基づいて、増幅された所定の標識付き核酸コピーの量の確率分布をモデル化し、
(ii) 標識付き核酸コピーの決定量が、前記標識付き核酸コピーのモデル化された確率分布に従って起こる尤度を決定し、
(iii) 推定増幅効率またはサンプル中の核酸コピーの期待数を変化させることにより、ステップ(ii)の尤度を最大化し、
(iv) 前記最大化モデルに従ってまたは前記最大化モデルにおける推定増幅効率に従って、サンプル中の核酸コピー数の推定値を提供する
各ステップを含む方法を提供する。

【0007】

本発明は更に、本発明方法のステップe)に従って、標識が取り付けられた核酸の増幅後の決定量に基づいたサンプル中の核酸コピーの量の推定値を計算するためのコンピューター読み取り可能装置を含むコンピュータプログラム製品を提供する。

【0008】

本発明の全ての実施形態および態様は組み合わせることができ、そしてあらゆる好ましい開示された実施形態は、本発明の全ての態様に対しても同様に参照する。すなわち、該方法の好ましい実施形態は、コンピュータプログラム製品の好ましい実施形態と解釈することもでき、またその逆も同様である。
〔発明の詳細な説明〕

【0009】

本発明は、上記要約の項目と特許請求の範囲において記載される。数種類の表現法が用いられるが、参照および実施例による概要の簡易性のために、実施例において用いられる多様な称号に関連して説明される。

【0010】

「核酸コピー数」および「核酸のコピーの数」（“Ｋ”）は、核酸コピーと称される共有の同一配列の核酸分子の量、存在度または濃度を指す。本明細書中で用いる場合、用語「コピー数」は、増幅前の、即ち元のサンプル中の量、存在度または濃度を示すために用いられる。これは、本発明の方法により決定することができる量またはコピー数である。そのようなサンプルは異なる量および／または配列の異なる核酸分子を有してもよい。この用語は異なる種の核酸に対しても呼称される。同一種の核酸分子（標識なし）は全て「プール」と称される。核酸分子の異なる「プール」は、異なる種の異なる核酸分子であり、各々の核酸種はそれ自身のバーチャルプールを形成する。「プール」は物理的分離を必要とせず、全ての核酸プールまたは種が同一サンプル中にあってもよい。

【0011】

「可変的標識」とは、識別することができる異なる性質または種の標識を言う。一例はランダムヌクレオチド標識、例えばランダムモノヌクレオチド（すなわちＡ、Ｔ（Ｕ）、Ｇ、Ｃの混合物）、ランダムダイマー（すなわちＡＡ、ＡＴ、ＡＧ等とＣＣの混合物）、ランダムトリマー、テトラマー、ペンタマー、ヘキサマー、セプタマー、オクタマー、またはより長鎖のヌクレオチドである。そのような標識混合物は様々な異なる標識を持ち、各々特異的な標識はおそらく１回または複数回表示されるだろう。ランダムヌクレオチドの場合、標識混合物の構造は、標識混合物プールのサイズに依存して、所定の標識の多重コピーに対して１つも存在しないかもしれない。標識は核酸分子に付着（結合）される。前記付着は非効率性に見舞われることがあり、これは全ての核酸分子が標識されるわけではないことを意味する。標識は指標「ｂ」により参照される。

【0012】

「核酸コピーに結合する標識の確率」（“ｐ_b”）は、核酸分子に標識を結合する確率として前記非効率性を数量化する。この数量化は、通常、本発明にかかる方法における推定である。

【0013】

「核酸複製手順」は、核酸の複製物を構築する手順である。そのような手法は、ＰＣＲのように、プライマーをアニーリングするステップと、前記プライマーを鋳型依存形式で伸長するステップを含んでよい。そのような手法は、非効率性に見舞われることがあり、すなわち全部の出発核酸分子が複製されるわけではないこと（すなわち「複製効率」が＜１であること、ここで１は、複製される各核酸分子の理論上の理想例であり、そして０は核酸分子が全く複製されないこと）を意味する。

【0014】

「推定増幅効率」（“ｐ_d”）は、前記増幅効率の推定値であり、本発明の計算手順において得られる。該推定値は、固有の増幅バイアス要因のために核酸種ごとに変動する、起こり得る増幅効率である（主要な異なる核酸配列）。本発明の推定モデルは、核酸種間の前記異なる増幅効率を反映する。

【0015】

「増幅された標識付き核酸コピーの量」（“ｎ_b（下付きｂ）”）は、増幅後の標識を含む核酸コピーの定量可能な量または濃度を指す。標識の種類は、変数および式の記述において下付き文字“ｂ”により示される。この表記は核酸コピーを指し、従って、同一のコピーを生成した１つの核酸種に関する。特に、その重要性を強調するために、ある種の核酸分子（または「プール」）には、それらに異なる標識が結合されていてよい。実際、標識結合効率を決定するためおよび次に本発明の方法の中での増幅効率をモニタリングするために、異なる標識を取り付けることが本発明方法において重要である。それらの異なる標識は、サンプル中のあらゆる種類の核酸に無差別に取り付けてもよい。「量」に関しての複数形は、標識を有する異なる増幅核酸コピーを測定できることを示す。上述した通り、「核酸コピー」は、同一の核酸種を指す。その違いは、個々の分子間で異なる標識によるものであるが、上述したような同じ標識種の多重標識も、核酸分子に取り付けることができる。これは、「各々の量が異なる標識を有する核酸コピーについて決定される」というように表現される。所与の核酸種の量（また、サンプル中のコピー数）を決定する推定手順によれば、“ｐ_b”（下付きｂ）”は、前記核酸コピー（前記種の核酸分子または「プール」）上の標識の量も指す。従って、一緒に研究される「前記標識の量」という表現は、標識を示す「ｂ」と一緒にも使用される（“ｎ_b（下付きｂ）”）。プール中の１つの標識付き増幅核酸コピーの量は、「Ｎ」と称される。それは標識「ｂ」（下付き）の個別量の総和（ｎ_b）である。

【0016】

「サンプル中の核酸コピー数の推定値」（「Ｋ」または「Ｋ_est」、ここで「est」は「estimated（推定）」の略であり、推定される性質を強調する）は、本発明方法によるサンプル中の核酸コピー数の推定値に関する。推定値はサンプル中の実際の自然量とは異なることがあるが、一方でそのような差または誤差は、従来技術の方法と比較すると大幅に減少される。本発明の方法の更なる有益性は、例えば、低い標識多様性および標識非効率性により導入される、標識バイアス（偏り）が大幅に減少されることである。「サンプル中の核酸コピーの期待数」という表現は、本発明方法が、各反復が新たな期待値を生成する期待値の精緻化（refinement）であり、最終的な「期待」値は推定値と称されることを示す。Ｋ_iまたはＫ_i+1 は、増幅前の標識付き核酸コピー数、即ち、標識付着後の核酸コピー数は、式中小文字の「ｋ」により示され、または標識種を表す「ｂ」を付けた「ｋ_b」により示される。この用語は通常、サンプル中の同一の塩基核酸配列（同一「プール」）を有するが、付着された標識が異なる核酸を指し示す。換言すれば、Ｋ_estは、計算に用いられる全ての標識に関して合計された１標識あたりの全核酸推定値〔ｋ_b〕の合計である。〔ｋ_b〕は次に標識確率〔ｐ_b〕に従った標識化の結果である（例えば、式(5) と(2)も参照のこと。Ｋ_est とｐ_b は式(3)を経由して相互に関連付けられる； Ｋ_est は式(3)への適当な解である）。

【0017】

「ステップd)での検出された異なる標識の数」（「ＤＵ」）は、異なる性質、例えば異なるバーコードまたは配列を有する検出された別個の標識の数を指す。従って、それは識別の数または種の型の数を指すが、各々の特異的な標識種の標識の量は意味しない。

【0018】

「増幅された標識付き核酸コピーの推定複製率」は、非効率性のために、増幅サイクルあたり各核酸分子を複製していない増幅反応におけるモデル化された複製率を指す。確率関数は結果として１サイクルあたりの統計的複製率に従うだろう。増幅（または複製）サイクル数は「ｃ」として示される。

【0019】

「標識付き核酸コピーのモデル化確率分布」（「ρ(ｎ_b)」は、標識付き核酸コピー（「ｎ_b」）に関する増幅過程の結果の上述したモデル化を指す。増幅後にある一定量に達する確率が計算される。様々な量についての集合確率が確率分布の中で統合される（量の数値に関して）。そのような確率分布関数は、ガウス分布形または他の形状をとることができ、標識付き核酸コピーが最大値のところにその量を有するであろう最大確率に近似した最大値を有することができる。分布関数の形状は、全体的最大と１以上の極大を有しうる。前記確率分布の生成は、推定増幅効率またはサンプル中の核酸コピーの期待数（それら２つは互換性である）に依存し、そしてそのモデルは、それら２つのパラメータのうちの１つの期待値を変更し精密化することにより調整することができ、それらのパラメータは、調査中の標識付き核酸コピー全体に渡って確率分布を実測量と最良適合（best fit）させるために尤度または確率を最大化する際に用いられる（サンプル中の全ての核酸である必要はなく全ての標識生成物である必要もない）。

【0020】

「増幅核酸コピーの決定量」（「Ｎ」）は、研究中の全ての核酸コピーに渡る総量を指す。「コピー」における複数形は、別個のコピーの数量的複数を指すことに注意。

【0021】

この概要に従って、本発明方法は、次のように、角カッコ内の指定要素の後に丸カッコ内の指定変数インジケーター（指標）と一緒に表記することもできる：
次のステップを含む増幅後の〔核酸コピー数〕（Ｋ）を推定する方法であって、
a) 決定すべき〔コピー数〕（Ｋ）の核酸を含むサンプルを提供し；
b) 可変的標識を前記核酸に付着させ；
c) 前記標識付き核酸を核酸複製手順を用いて増幅させ；
d) 〔増幅された標識付き核酸コピーの量〕（ｎ_b）を決定し、その各々の量が、〔異なる〔標識〕（ｂ）を有する核酸コピー〕について決定され（異なるｂについてのｎ_b）；
e) 以下により決定された〔ステップd)の量〕（ｎ_b）に基づいてステップa)の〔サンプル中の核酸コピー数の推定値〕（Ｋ_est）を提供し、ここでＫ_estは
サンプル中の少なくとも２、好ましくは少なくとも４、少なくとも10以上の異なる核酸種について、〔１つの核酸に結合する標識の確率〕（ｐ_b）を概算することにより決定され；前記ｐ_bは増幅後の異なる核酸種について〔前記標識の量〕（ｎ_b）を平均化し；そして反復的または段階的に
（Ａ）
〔ステップd)において検出された異なる標識の数〕（ＤＵ）、〔核酸コピーに結合する標識の確率〕（ｐ_b）に従った核酸に結合する異なる標識の数の期待値および〔サンプル中の核酸コピーの期待数〕（Ｋ_est）または〔推定増幅効率〕（ｐ_d）の事前反復もしくは既定値に基づいて、〔サンプル中の核酸コピー数の推定または期待数〕（Ｋ_est）を精密化（改善）し；または
（Ｂ）
(i) 決定された〔ステップd)の量〕（ｎ_b）、〔１つの核酸に結合する標識の確率〕（ｐ_b）、〔推定増幅効率〕（ｐ_d）または〔サンプル中の核酸コピーの期待数〕（Ｋ_est）に基づいて、そして前記〔推定増幅効率〕（ｐ_d）に依存する複製サイクルあたりの核酸複製手順における増幅された〔標識付き核酸コピー〕（「プール」）の推定複製率に基づいて、〔所定の標識を有する増幅された核酸コピーの量の確率分布〕（ｐ(ｎ_b)）をモデル化することにより概算し、
(ii) 〔標識付き核酸コピーの決定量〕（ｎ_b）が前記〔標識付き核酸コピーのモデル化確率分布〕（ｐ(ｎ_b)）に従って生じる尤度を決定し、
(iii) 〔推定増幅効率〕（ｐ_d）または〔サンプル中の核酸コピーの期待数〕（Ｋ_est）を変化させることにより、ステップ(ii)の尤度を最大化し、
(iv) 前記最大化モデルに従ってまたは前記最大化モデルにおける〔推定増幅効率〕（ｐ_d）に従って、〔サンプル中の核酸コピー数の推定値〕（Ｋ_est）を提供する
ことにより、決定される
ステップを含む方法。

【0022】

本発明方法は、決定すべきコピー数の核酸を含むサンプルをベースにする。該サンプルは、様々な異なる核酸種を含んでよく、各核酸種が１以上の分子または「コピー」により表される。例えば、サンプルは1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 20, 30, 50, 100, 150, 200 以上、またはそれらの数値の間の任意範囲の数の異なる核酸種を含んでよい。少なくとも１以上の種のコピー数は、1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 20, 30, 50, 100, 150, 200, 1000, 10000, 20000 またはそれ以上の分子またはそれらの数値の間の任意範囲のコピー数であることができる。

【0023】

核酸は全長の天然核酸分子、例えばmRNA、ミクロRNA、rRNA、tRNAまたはゲノム核酸もしくはベクター核酸、例えばゲノム（好ましくは短鎖のもの、例えばウイルスもしくは細菌ゲノム、またはそれらの分子成分）、プラスミドもしくはトランスポゾンのようなベクター、または人工DNA構成物であることができる。一般に、核酸はRNAまたはDNAであるか、あるいは増幅可能である任意の他の核酸型であってよい。

【0024】

好ましくは、核酸は、大型の核酸由来のフラグメント、例えば次世代シークエンシング（NGS）またはRNA-Seqまたはショットガンシークエンシングにおいて生産されるフラグメントである。フラグメント化は、せん断のような物理的手段により、または制限酵素切断などの化学的または酵素的手段によることができる。フラグメント化は、無作為（ランダム）、すなわちランダムフラグメントを生じるように切断の部位特異性が全く存在しなくてもよく、あるいはそれは部位特異的、例えば酵素制限消化のように、選択された核酸パターンに対して特異的であってよい。ランダムフラグメント化は、多数の異なる核酸種を生じるが、それは計算法によりアシストされる最新のシークエンシング法により処理することができる。

【0025】

好ましくは、本発明方法により分析しようとするサンプルの核酸（例えば上述したようなフラグメントまたは他の核酸分子）は、10～10,000ヌクレオチド（nt）、好ましくは15～8000、または20～5000、または50～4000、または80～3000 ntの平均長さを有する。

【0026】

ステップb)において、可変標識が前記核酸に付着される。付着は、既知の化学的または酵素的反応、例えばライゲーション反応または結合反応により促進することができる。付着は、ステップd)の検出まで前記標識が操作の間ずっと核酸に付着されたまま保持されるのを促進する。核酸を標識するのに適当である任意の標識を選択することができる。標識は、ステップc)の増幅反応において増幅生成物にも結合することが必要である。従って、最も便利な標識は、バーコードまたはUMIのような核酸タグである。しかしながら、追加の労力により、タンパク質やペプチド、ペプチド性ヌクレオチド、抗体、レセプター、抗原、認識分子（ビオチンとアビジン、蛍光標識などの結合パートナーにより識別することができる）、量子ドット等のような他の標識も使用することができる。好ましいのは、高変動性を可能にする標識である。同じこの理由により、ヌクレオチド標識が好ましいが、それらは唯一の可能性ではない。

【0027】

好ましくは、標識は増幅する能力があり、特に、国際公開WO2017/051387 A1パンフレットに開示されるように増幅させることができる巨大分子である。

【0028】

ヌクレオチド標識の場合、それらはRNAもしくはDNAまたは任意の他のヌクレオチド型であることができる。それは核酸分子と同じ型または他の型であってよい。核酸標識は、認識要素として標識分子の固有の特性を利用する。好ましくは、前記認識要素がヌクレオチド配列である。ヌクレオチドは、Ａ、Ｔ（Ｕ）、Ｇ、Ｃまたはそれらの任意組み合わせ、このましくは４つのヌクレオチド型全て（ＵはRNAにおいて優先的に使用され、ＴはDNAにおいて優先的に使用され；両者は相補的Ａにより認識される）から選択され得る。標識は１，２，３または４種の異なるヌクレオチド型を含んでよい。唯一のヌクレオチド型を有する場合、認識要素は必然的に標識のサイズまたは長さである。複数のヌクレオチド型を含むものは、長さが６ヌクレオチドである４種の異なるヌクレオチド型を有する標識について、多数の順列（並び替え）、例えば４⁶（４exp6）が可能である。順列の数は、ヌクレオチド型の数の後に（長さnt）の指数を付ける定法、例えば4³, 4⁴, 4⁵, 4⁶, 4⁷, 4⁸ 等、または３種のヌクレオチド型の場合には3³, 3⁴, 3⁵, 3⁶, 3⁷, 3⁸ 等に従う。異なる標識の適当数は、目的のサンプルに予想される核酸の複雑性と多様性に従って、実施者により決定することができる。例えば、2, 3, 4, 5, 6, 7, 8, 9, 10, 15, 20, 25, 30, 35, 40, 50, 60, 70, 80, 100, 120, 140, 160, 180, 200, 300, 400, 500, 750, 1000, 2000, 5000, 10000もしくはそれ以上、またはそれらの数値の間の任意範囲の異なる標識数が使用される。幾つかの方法では、特にメソッド（Ａ）では、好ましくは異なる標識の数は、それが推定核酸コピー数の数よりも大きくなるように選択される。このことはメソッド（B)でも可能であるが、必然ではない。メソッド（Ｂ）は、より少ない数の異なる標識を取り扱うのにより有効である。

【0029】

標識は必ずしも分子として提供されなくてもよいが、多様なコピー数の標識を有することも可能である。本発明方法は、任意の標識バイアス、即ち可変的標識の不均一濃度および核酸分子への標識の不均一な付着が原因の任意のバイアス（偏り）を相殺する能力を有する。好ましくは、個々の標識種は、1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 50, 100, 1000, 10000 またはそれ以上、例えば好ましくは100万（M）まで、100,000 までまたは10,000までのコピー数を有する。好ましくは、全標識の総数は、調査中の各核酸が標識されるためには全ての種の全核酸の総数に等しいかまたはそれを超えるべきである。全標識の総コピー数は合わせると少なくとも1000 または10000またはそれ以上、例えば好ましくは10000百万 (M)まで、5000 Mまで、1000 Mまで、または500 Mまでであってよい。

【0030】

好ましい実施形態では、標識は、標識配列中のエラー（誤り）の検出と補正を可能にするものが使用される。そのようなエラー補正標識はそれ自体当業界で既知であり、Krishnan 他、Electronics Letters, 2011, 47, 236-237に記載されている。それらの標識は、標識配列のシークエンシングの間の配列決定の誤り－所定のシークエンサーの許容誤差を検出するのに有用である。例えば、配列ACGTとTGCAを有する２つの標識を考慮する場合、１個のヌクレオチドが別のものに置き換えられているそれらの標識中の１つのエラーを補正することができる。何故なら、それらの標識は４つの位置全てが異なるためである。例えば、もし第一の標識の読み出しがACCTに変更され、そして第二の読み出しがTCCAに変更されるならば（両方とも配列決定の誤りによる）、その正しくない（誤った）標識は、それらの正しい形から１の距離を有するが、別の標識の正しい形からは３の距離を有する。よって、正しい標識からのそれらの距離に基づいて、ACCTはACGTと正しく割り当てられ、TCCAはTGCAに割り当てられるだろう。他方で、もし第一と第二の標識が両方ともACCAに変更されたとしたら、これはその両者の正しい標識に対して２の距離を有するだろう。よって、このエラーを補正することができず、１より大きいエラーが発生したと推測するだけである。結果として、標識ACGTとTGCAは１エラー補正および２エラー検出用である。この例は、置換エラーの修正を参照する（Krishnan他、前掲; およびBystrykh他、PLOS ONE, Public Library of Science, 2012, 7, 1-8; 全ての参考文献は参照により本明細書中に組み込まれる）。全ての異なる標識間に大きなヌクレオチド差異を有する長鎖の標識は、より多数のエラーの補正と検出を可能にする。置換、挿入および欠失を補正する標識を作製することも可能である（Buschmann他、BMC Bioinformatics, 2013, 14, 272; Hawkins他. Proceedings of the National Academy of Sciences, National Academy of Sciences, 2018, 115(27):E6217-E6226; 全ての参考文献は参照により本明細書中に組み込まれる）。好ましくは、そのようなエラー補正標識は、本発明方法のステップb)において使用される。特に好ましくは、エラー補正標識が置換－、挿入－および欠失－補正標識である。特に好ましくは、それと組み合わせてまたは代わりとして、エラー補正標識が１，２，３またはそれ以上のエラーを補正するのに適当である。そのような標識は、２ｎ＋１だけの配列の差異を有することができ、ここでｎは補正可能であるエラーの数である。好ましくは、配列の差異が２ｎ＋２である。式２ｎ＋１に比較した前記追加の差異は、上記の例に記載のようにエラー補正量の上に別のエラー検出の層を追加するだろう。

【0031】

好ましくは、ステップd)において、エラー補正は、正しい標識（配列決定の誤りのない）を割り当てるのに用いられる。純エラー検出（補正なし）は、更なる分析からそのような標識を取り除くために用いられるだろう。

【0032】

ステップc)は、PCRのような核酸複製方法により標識核酸を増幅させることを含む。核酸増幅方法は、一般に、所定の鋳型で標識された核酸について、（更なる）コピーが生成されるというパターンに従う。そのような方法は、プライマー結合と前記プライマーの伸長を含むことができる。プライマー結合のステップは、リンカーまたはアダプター分子をプライマー結合領域に付着させることを必要としうる。これは、先端または末端から出発して任意の核酸配列を増幅することができるという利点を有する。別のあまり好ましくない配列非依存性プライミングは、ランダムプライミングである。発明の背景の項目において記載した通り、核酸増幅方法は非効率性、および前記非効率性の配列バイアス（偏り）にも悩まされる傾向がある。このバイアスは一般には知られていない。本発明は、本質的にそのような未知のバイアスと非効率性を埋め合わせる回避方法を提供する。

【0033】

増幅の様式、すなわち一般的な場合での複製は、本発明の計算ステップにおいてモデル化され考察される。増幅の様式が異なり得る場合、例えば増幅サイクルあたり三重複製または他のｘ倍複製である場合、これは本発明方法において同様に上手く検討することができる。PCRのような標準増幅法は、サイクルあたり（二重）複製であるため、本発明は複製に関して記載されている。もちろん、本発明は代わりにサイクルあたり別のｘ倍増幅とも読みとれる。もちろん、複数サイクルの増幅が通常使用される。従って、複製方法は何倍もの増幅を生じうる。好ましくは、2, 3, 4, 5, 6, 7, 8 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30 またはそれ以上の増幅サイクルが使用される。

【0034】

ステップd)は、増幅された標識付き核酸コピーの量を決定することを含み、各々の量は異なる標識を有する核酸コピー毎に決定される。量の決定は、一般に、混合物中の核酸量を推測するための標識付き核酸の同定または分離を含む。恐らくおよび通常は分子ごとに異なる標識は、核酸種のアイデンティティーと一緒に同定される。代表的方法としては、核酸分子の配列と共に、それのアイデンティティー（ID）も提供し、次いでカウントすることができる配列決定法が挙げられる。従って、このステップは、増幅後の核酸コピーの量（標識と一緒に）を提供する。ここでは、ステップe)で実施される増幅以前の核酸の量または数に到達することが目標である。

【0035】

ステップe)によれば、本発明方法は、（Ａ）および（Ｂ）と名付けられた２つの選択可能な具体的メソッドを提供し、これらはそれぞれ実施例2.1.2および2.2により更に裏付けられる。一般的には、メソッド（Ａ）および（Ｂ）の両方について、ステップa)のサンプル中の核酸コピー数の推定量は、サンプル中の少なくとも２、好ましくは少なくとも４、少なくとも10またはそれ以上の異なる核酸種について：増幅後の異なる核酸種ごとに前記標識の量を平均化することによって、１つの核酸に結合する標識の確率を概算することにより、ステップd)の決定量に基づいて提供される。前記標識の量を平均化するステップは、再度、本発明方法において用いられる可変的標識について異なる核酸種（またはフラグメント）に渡り平均化された標識の濃度のコピー数を参照する。従って、特異的標識が１つの核酸分子にどれくらいの頻度で付着されたか―それは平均化されるため核酸種に依存しない―を決定することができる。異なる核酸種への標識結合が異なり得ると仮定すれば、これは確率的影響であり統計的変動性に従うため、標識の結合不均一性と結合反応の不均一性により更に影響を受けるだろう。平均値を得るためには、少なくとも２つの核酸種、好ましくは２より多く、例えば3, 4, 5, 6, 7, 8, 9, 10, 20, 30, 40, 50もしくはそれ以上の異なる核酸種またはそれらの数の間の任意範囲に対して平均化が実施される。次いでこの平均値は、標識結合確率（「核酸コピーに結合する標識の確率」とも称される）の決定においてまたは核酸に結合する異なる標識の数の期待値の決定において用いることができる。そのような確率または期待値は、前記平均および特定の核酸種（すなわち平均化されたメンバー）についての前記平均からの偏差を考慮に入れて、標準的な統計学的方法に従って算出することができる。これは、標識の数が、標識結合反応速度論および標識混合物組成のどんな不均質性（恐らく他のものと比較して或る種の標識が過剰発現されている）にも依存することを意味する。増幅効率がいったん既知になれば（所定の核酸種について）、増幅後に決定された量を使ってサンプル中の核酸の元の数を直接算出することが可能である（実施例の項目の式（27）を参照のこと）。増幅効率が直接導出できない場合、前記増幅効率を反復して概算することが可能である。増幅効率は核酸数に直接関連するので、前記数は、増幅効率の代わりにまたはそれに加えて、同等に計算に送出する（フィードする）ことができる。事実、好ましい実施形態では、推定増幅効率またはサンプル中の核酸コピーの期待数は互換的に使用される。それらは式（27）に従って変換することができる。よって、本発明方法が「推定増幅効率」という表現を使用する時は常に、前記表現は「サンプル中の核酸コピーの期待数」も指し、そしてその逆も同様である。

【0036】

好ましくは、前記核酸複製手順において前記標識付き核酸を増幅する効率は100％ではない。それは異なる種の核酸コピーごとに異なり得る。それは本発明方法が本質的に効力入れて補償する難しさである。通常の増幅効率は60％～99.99％、好ましくは65％～99％、例えば70％～95％、最も一般的には約80％である。核酸コピーの異なる種間の代表的平均変動は、最大±20％または最大±10％である。

【0037】

本発明は、反覆プロセスにおけるサンプル中の核酸コピー数の推定を提供する。これは、核酸コピー数の一次的な推定または概算が段階的に改善されることを意味し、通常、核酸コピー数の最終推定に達するまで繰り返し改善される。段階的な改善の代わりに、様々な核酸推定値を計算する「総当たり（brute-force）」アプローチも可能であり、その場合、様々に変動する核酸推定値が（可能な全範囲または合理的範囲、例えば１からステップd)で決定された全ての併合した核酸種量の総量（“Ｎ”）までをカバーする範囲において）、間隔を空けてまたは段階的に計算され、そしてこれらすべての推定値の最善モデルが選択される。反復プロセスは、実施者がその改善に満足した場合、実施者が停止させることができ、それは通常、収束後または最適適合（best fit）の選択により十分な数の期待値が推定された場合である。収束は、例えば、実施者が満足した時、例えば反復によって核酸コピー数の推定値が変更されないか、あるいは５％以下、例えば３％以下、または１％以下、または0.5％以下、または0.3％以下だけ変更される時、実施者により設定される。本発明方法で使用される他のパラメータは、特定の増幅モデルに適合するための標識付き核酸コピー（増幅後）の決定量など、反復を停止する決定にも同様に使用できる。次いで、前記モデルを使用して、サンプル中の核酸コピー数を算出することができる。本発明は、サンプル中の核酸コピーの最大推定数に向けた反復的な改善（精緻化）を含む。多くの最大化手順が当技術分野で知られており、本発明に従って使用することができる。考えられる反復手順の例は、０から１の増幅効率区間を等間隔の部分区間に分割する。これらの部分区間から、選択された部分区間の境界点での尤度が未選択の部分区間の境界点での尤度よりも高くなるように、１つ以上の部分区間が選択される。選択された各部分区間について、区間内で最尤度を有する増幅効率を見つける二分探索[4]が実行される。全体の最大値は、二分探索結果の最大値である。部分区間の数は１からＮまでであることができ、二分探索は１からＮ個の部分区間で実行できる。もちろん、他の多くの最大化方法も可能である。同じく、サンプル内の核酸コピー数の期待数（予想数）を、推定増幅効率の代わりに使用することができ、これは互換的に使用できる。例では、Ｋ_estは式(3)を（ほぼ）満たす値である。ここで、(3)の右辺は(2)で定義される。(3)のような非線形最適化または非線形プログラミングが存在する。一例は式(4)および(5)に従った不動点反復である。

【0038】

メソッド（Ａ）と称される本発明の特定の実施形態は、サンプル中の核酸コピーの推定数または期待数（すなわち、反覆後に最終推定値が決定されるまでの中間）を精緻化することを含む。前記精緻化（refining）は、ステップd)で検出された個別の標的の数、および１核酸コピーに結合する標識の確率とサンプル中の核酸コピーの期待数の既定値または事前反復または推定増幅効率（上記のように、これら２つの「または（or）」変数は互換的に使用できる）に従った１核酸に結合する異なる標識の数の期待値に基づくまたはそれを使用する。ステップd)で検出された個別の標識の数は実測値であり、更に説明する必要はない（他の場合は、上記の用語定義の項目とステップd)を参照されたい）。１核酸コピーに結合する標識の確率に従った１核酸に結合する異なる標識の数の期待値、またはより良い参照のための可変的用語：[１核酸コピーに結合する標識の確率]（ｐ_d）または「標識結合確率」で記述されるものは、本質的に標識効率の期待値を指し、これは増幅後の異なる核酸種についての前記標識の量の平均に基づいて上記のように計算することができる。１核酸コピーに結合する標識の確率は、それの相関によって、すなわち前記確率が前記期待値と直接の相関があるため、１核酸に結合する異なる標識の数の期待値を算出するために用いることができる。統計および確率演算に熟練した数学者は、そのような計算を行うことができる。一例は、本明細書中の実施例の項目に式(2)で与えられる。また、ステップe)の概論に関する上記項目を参照のこと。最終的なパラメータは、サンプル内の核酸コピー数の期待数の既定値（プリセット値）もしくは推定増幅効率、または組み合わせモデルではその両者である。最初の反復は通常は既定値から出発する。通常、その後の反復は、直前の反復の値から開始することができる。通常、サンプルを提供する実施者は、推定値を見事に推測し、それを本発明の方法で利用できる。その他の適切な開始値は、一般的または通常の増幅効率をベースにする。例えば、サンプル中の核酸コピーの期待数の初回反復の既定値、またはステップd)で検出された異なる標識の数の値の範囲内のまたは１から増幅核酸コピーの決定量までの範囲の整数による推定増幅効率を選択することが可能である。この実施形態はまた、より好ましく参照するために、実施例で使用されるカッコ内のパラメータ表現および変数と共に次のように記述される：[サンプル内の核酸コピーの期待数]（Ｋ_est）または[推定増幅効率]（ｐ_d）の初回反復の既定値は、[ステップd)で検出された異なる標識の数]（ＤＵ）の値から、または１から[増幅された核酸コピーの決定量]（Ｎ）までの範囲の整数によって選択される。このメソッドは反復的改善を提供するため、開始値または既定値はそれほど重要でない。Ｋ_est値は継続的に改善されるだろう。遠く離れた開始既定値は、より適した開始値または既定値で開始するプロセスと同じ最終品質に到達するために、さらに反復する必要があるだろう。もう１つの選択肢（Ａ）のメソッドは、実施例2.1.2およびその中に示される式において更に例示され、そのメソッドは本明細書に記載の一般的発明に従って使用することができる。

【0039】

本発明者らは、サンプル中の核酸コピー数を計算するための更により改善された方法（メソッド）も提供している。このメソッドも反復的であり、異なる核酸種についての前記標識の量の平均、または核酸に結合する標識の確率を使用する。

【0040】

メソッド（Ｂ）は、(i)から(iv)と称する４つのステップを本質的に含む。このメソッドは、通常、(i)所定の標識付き増幅核酸コピーの量の確率分布をモデル化することを含む。個々の核酸種の量の総和、すなわち、研究中の全ての核酸種（「プール」）についての標識付き増幅核酸コピーの量を代替的にまたは付加的に使用することも可能である。上記の用語定義の項目を参照のこと。簡単に言えば、増幅プロセスのモデルが作製され、このプロセスは、増幅の機構を考慮に入れ、次にサイクルあたりの理論上の増幅率（通常は複製）、サイクル数および増幅効率に依存し、そして湿式化学反応の増幅効率が不明であると仮定してそれらの挙動をモデル化する。このメソッドは、代わりに、ステップd)の核酸（標識付き）の決定量（既知の決定値）、核酸に結合する標識の確率（上記で決定された近似値；これはメソッド(Ａ)と共通のステップである）、推定増幅効率またはサンプル中の核酸コピーの期待数（上記のように互換的に使用され、それらの値は再度反復的に改善される）を使用し、そして前記推定増幅効率に依存する複製サイクルごとの核酸複製手順における標識付き増幅核酸コピーの推定複製率に基づく（すなわち増幅プロセスのモデル）。そのようなモデルは、一般に当技術分野で知られている。増幅方法は広く研究されており、（非効率的）複製またはその他のｘ倍濃縮工程における既知のサイクル数の後の統計的増幅率が研究されている。参考文献[1], [2], [3]および実施例の項目2.2を参照のこと。好ましくは、標識付き増幅核酸コピーの推定複製率は、ガウス分布、負の二項分布、ガンマ分布、ディラックのデルタ分布、もしくはガルトン・ワトソン分布の複製確率関数、またはそれらの混合分布でモデル化される。好ましくは、所定の標識付き増幅核酸コピーの量の確率分布をモデル化する際に、確率的分岐モデルが使用される。このようなモデルは、ガルトン・ワトソン（ＧＷ）モデルなどの複製増幅プロセスの性質を考慮する。計算時間要件が低いモデルは、例えば負の二項分布または単峰型分布または多峰型（マルチモーダル）分布である。また、混合分布を使用して、より高位の数の標識付き核酸コピー（「ｋ」）値についての確率関数を提供するといった複雑なタスクのために、低計算時間要件を有するモデルを使用するなど、コンピュータの計算時間を短縮することができる。サンプルに存在する前記「ｋ」値は未知である場合もあるが、「ｋ」値は式(25)または(7)に従う変数として、ｐ（ｋ|Ｋ,ｐ_b）により重み付けされたｐ（ｎ_b|ｋ）の総和として近似することができる。後者は、式(26)または(8)で与えられる。確率分布の変数を加算することにより、この変数が確率分布から取り除かれる。これは通常、未知の変数に行われる。式(25)によると、ｐ（ｎ_b|Ｋ）は、ｐ（ｎ_b|ｋ）と重みｐ（ｋ|Ｋ,ｐ_b）の混合物である。方程式(25)と(26)は、方程式(8)と(7)の特殊ケースである。例えば、ｋについての高位数とは５以上、好ましくは10以上、例えば20、30、40、50、60以上、またはこれらの値の間の任意範囲のｋである。そのような確率関数は、様々な目的の核酸種について計算される。従って、多様な確率分布が取得される。

【0041】

メソッド（Ｂ）のステップ(ii)は、ステップ(i)の標識付き核酸コピーのモデル化確率分布に従って起こる、標識付き核酸コピーの決定量（ステップd)におけるような）の尤度を決定することを含む。再び、「推定増幅効率」と「核酸コピーの期待数」は互換的に使用できることが強調される。従って、変動を必要とするパラメータはただ１つである。数値は、例えばランダムにまたは決定的に、例えは反復的なステップを踏むことなどの更なる情報を利用することによってまたは確率分布関数の情報を利用することによって、例えば数値を予想される方向に移動させることにより、例えば決定値が確率分布関数の最大値の方向に移動するように曲線の山（上り坂）の上方に移動させることなどにより、変動される。最大化プロセスのための多くのメソッドが当技術分野で知られており、本発明に使用することができる。

【0042】

最後に、ステップ(iv)において、サンプル中の核酸コピー数の推定値が、ステップ(iii)の前記最大化モデルに従って提供される。最大化は、例えば実施者がサンプル中の核酸コピー数の改善に満足した時に、またはステップe)の導入部において上述した通り、収束に到達した時に、反覆プロセスに至る。これでステップe)が終結する。

【0043】

本発明は、ステップe)の結果を、プリンタ出力、スクリーンなどの可読媒体上に提示すること、またはハードドライブ、フラッシュメモリなどのコンピュータ記憶装置などのデータキャリアに書き込むことを更に含みうる。

【0044】

一般に、本発明は、実施例およびそこに示される式で更に説明される。これらの式を個別に使用して、本発明を説明し、本発明方法を更に特定することができる。例えば、標識が１核酸に結合する確率（「ｐ_b」）は、式(12)、式(45)または式(46)に従って決定できる。推定増幅効率（「ｐ_d」）またはサンプル中の核酸コピーの期待数（「Ｋ」または「Ｋ_est」）を互換的に使用でき、そして式(27)に従って変換することができる。メソッドＡ）では、所定の標識付き増幅核酸コピーの量の確率分布（「ｐ(ｎ_b)」）は、式(3)に従って、好ましくは式(4)と(5)の両方によって決定することができる。もう１つの選択肢メソッドＢ）では、所定の標識付き増幅核酸コピーの量の確率分布（「ｐ(ｎ_b)」）は、式(25)に従って決定することができる。これらの式の利用はすべて、実施例と同様に本発明の好ましい実施形態であり、それらは互いに組み合わせることができる。

【0045】

本発明は、本発明方法を使用するコンピュータプログラム製品を更に提供し、例えば、コンピュータ上で前記メソッドおよびステップを実装または支援するための機械語を含む。コンピュータプログラム製品は、あらゆる種類の記憶装置上に提供できる。また、本発明の方法のステップの実装を支援するようにプログラムされたシステム、例えばコンピュータ装置も提供される。計算ステップは通常、オペレータの補助なしで実行される。入力および設定ステップは、プログラムまたはシステムによって、例えば必要に応じて、ランダムステップの繰り返し回数に対するオプション提案を提示することによって支援される。もちろん、プログラムまたはシステムは、オペレータからの追加の入力無しにデフォルトパラメータを使って実行されてもよい。特に、本発明は、ステップe)に従った標識付き核酸の増幅後の決定量に基づいて、サンプル中の核酸コピーの量の推定値を計算するためのまたは計算に適合させた、コンピュータ可読命令を含む、コンピュータプログラム製品を提供する。言い換えれば、本発明は、コンピュータによって実行されると、コンピュータに本発明のメソッドの少なくともステップe)を実行させる命令を含むコンピュータ可読媒体を提供する。もちろん、他のメソッドのステップのデータを上記に詳述した通りステップe)で使用可能である。

【0046】

コンピュータプログラム製品を備えるコンピュータ可読記憶装置は、コンピュータ上で本発明のメソッドを実装するように、またはコンピュータによって本発明のメソッドを支援するように適合されている。特にステップe)はコンピュータで実行される。通常の湿式化学であってもステップa), b), c)および／またはd)は、例えば自動または半自動シークエンスリーダー（配列解読装置）からのデータを制御および取得するために、コンピュータにより支援され得る。コンピュータプログラム製品または記憶装置は、シークエンサー、好ましくはコンピュータ構成部品を含むシークエンサーなど、サンプルから短いシークエンシングリードを取得するシーケンス生成部品を備えてもよい。例えば、コンピュータ可読媒体としては、磁気記憶装置（例えばハードディスク、フロッピー（登録商標）ディスク、磁気ストリップなど）、光ディスク（例えばコンパクトディスク(CD)、デジタル多用途ディスク(DVD)など）、スマートカード、およびフラッシュメモリ装置（カード、スティック、キードライブなど）が挙げられるが、それらに限定されない。ステップa) ～d)の実行のための適応は、コンピュータプログラム製品の一部であるか、または一部ではない。コンピュータプログラム製品が、決定量の入力情報を受理できれば十分である。標識付き増幅核酸コピーの各量は、ステップd)に従って、異なる標識を持つ核酸コピーに対して決定される。従って、コンピュータプログラム製品は、前記入力からステップe)を実行するように適合されている。また、増幅の性質（例えばPCRのような複製方法）および好ましくは増幅のサイクル数および潜在的に結合した標識が、コンピュータプログラム製品の入力情報として使用される。

【0047】

従って、本発明は、特にコンピュータが前記入力を受信した後、コンピュータでステップe)を実行することにも関する。

【0048】

コンピュータプログラム製品は、ステップe)の結果を、プリンタ出力、スクリーンなどの可読媒体上に提示する（書き込みを含む）または表示するように適合させる、またはハードドライブ、フラッシュメモリなどのような、例えは上述したような、コンピュータ記憶装置等のデータキャリア上に書き込まれるように適合させることが可能である。同様に、コンピュータ上でステップe)を実行する方法は、好ましくは、斯かる可読媒体上にステップe)の斯かる結果を提示するステップを更に含む。

【0049】

本発明は、本発明のこれらの実施形態に限定されることなく、以下の図面および実施例によって、さらに説明される。特定の実施形態では、本発明は、不均等分布した標識のセット（バーコード、UMIなど）または単純標識の計数（カウンティング）には小さすぎる標識のセットについて、PCR前（プレPCR）（または他の増幅）核酸数の精度が高い推定値を生成するメソッドを提供する。この目的のために、最初に増幅前の標識分布を推定し（図４）、その後、それをPCR前の核酸数の推定に使用する。第一に、本発明者らは、増幅後に観察された異なる標識の数を入力情報として使用するだけでなく、増幅前の標識分布も考慮に入れた、標識計数の改良を研究する。残りのメソッドは、入力としての増幅後標識度数を頼りにし、統計モデルを利用して増幅プロセスを説明する（図５）。例として、PCRプロセスのポアソンモデル、二項分布モデル、およびガルトン・ワトソンモデル用に作成された合成データに関して本発明メソッドを評価する。

【図面の簡単な説明】

【0050】

【図1】配列特異的なPCRバイアス。同じサイズの同一cDNAフラグメントのプール（cDNA1およびcDNA2と標識付け）をPCRにより増幅した。PCR後、それらのサイズは異なる。これは不正確な遺伝子とイソ型の定量化に至る可能性がある。

【0051】

【図2】UMIはPCRにより生成されたフラグメントコピーの識別を助ける。UMIはPCR前のcDNAフラグメントに連結される。UMIのセットが十分に大きいならば、同じUMIが１cDNAフラグメントの２コピーに連結することはないだろう。この場合、PCR後の個別のUMIの数は、PCR前のフラグメントコピー数と同じである。

【0052】

【図3】PCR前プールサイズの正確な推定は、均等分布したUMIの大型セットを必要とする。UMIセットが小さすぎるかまたは分布が均等でない場合、同じUMIが２つの異なるフラグメントコピーに連結しうる。PCR後に個別のUMIをカウントすることは、PCR前コピー数の過小評価につながる可能性がある。

【0053】

【図4】各プール中のPCR後UMI分布からUMIのPCR前分布が推定される。このためには、全フラグメントプールに渡りPCR後のUMIカウントの度数が平均化される。

【0054】

【図5】プール中のPCR前フラグメント数の推定。モデル依存性メソッドは、入力データとしてPCR後のUMIカウントを用いる。モデル非依存性メソッドは、入力データとして個別のUMIの数を使用する。どちらのタイプのメソッドも、入力データとして推定PCR前UMI分布を使用する。

【0055】

【図6】確率分布

【数1】

（実線）対、負の二項分布（破線）および正規分布（点線）による近似。各パネルはｋ＝1,....,5 についての分布を示す。左から右方向へのピークは、昇順でのｋに対応する。

【0056】

【図7】UMI確率ｐ_bの確率密度関数。各UMI位置においてα＝５を有するディリクレ分布から（ｐ_A, ｐ_c, ｐ_G, ｐ_T）が抽出される。UMI確率ｐ_bは、そのヌクレオチドの確率の積として与えられる。オレンジ色の垂線は、ｘ軸上の１／Ｂの所に引かれ、均等UMI分布

【数2】

の場合はｐ_bの値である。

【0057】

【図8】UMIカウント値の精度。二項分布モデルからのデータ。UMI分布のタイプ（均等／不均等）およびUMIの数Ｂがグラフの表題に示される。複製前フラグメントの推定数Ｋの相対誤差がｙ軸上に、真の数がｘ軸上に示される。

【0058】

【図9】度数補正したUMIカウントの精度。二項分布モデルからのデータ。UMI分布のタイプ（均等／不均等）およびUMIの数Ｂがグラフの表題に示される。複製前フラグメントの推定数Ｋの相対誤差がｙ軸上に、真の数がｘ軸上に示される。

【0059】

【図10】ポアソン－正規分布メソッドの精度。Ｂ＝64 UMIを用いたポアソンモデルからのデータ。UMI分布のタイプ（均等／不均等）および効率ｐ_dがグラフの表題に示される。複製前フラグメントの推定数Ｋの相対誤差がｙ軸に、真の数がｘ軸上に示される。

【0060】

【図11】ポアソン－正規分布メソッドの精度。Ｂ＝256 UMIを用いたポアソンモデルからのデータ。UMI分布のタイプ（均等／不均等）および効率ｐ_dがグラフの表題に示される。複製前フラグメントの推定数Ｋの相対誤差がｙ軸に、真の数がｘ軸上に示される。

【0061】

【図12】二項－正規分布メソッドの精度。Ｂ＝64 UMIを用いたポアソンモデルからのデータ。UMI分布のタイプ（均等／不均等）および効率ｐ_dがグラフの表題に示される。複製前フラグメントの推定数Ｋの相対誤差がｙ軸に、真の数がｘ軸上に示される。

【0062】

【図13】二項－正規分布メソッドの精度。Ｂ＝256 UMIを有する二項分布モデルからのデータ。UMI分布のタイプ（均等／不均等）および効率ｐ_dがグラフの表題に示される。複製前フラグメントの推定数Ｋの相対誤差がｙ軸に、真の数がｘ軸上に示される。

【0063】

【図14】ＧＷ多成分混合の精度。Ｂ＝64 UMIを用いたガルトン・ワトソン（Galton-Watson）モデルからのデータ。UMI分布のタイプ（均等／不均等）および効率ｐ_dがグラフの表題に示される。複製前フラグメントの推定数Ｋの相対誤差がｙ軸に、真の数がｘ軸上に示される。

【0064】

【図15】ＧＷ多成分混合の精度。Ｂ＝256 UMIを用いたガルトン・ワトソン（Galton-Watson）モデルからのデータ。UMI分布のタイプ（均等／不均等）および効率ｐ_dがグラフの表題に示される。複製前フラグメントの推定数Ｋの相対誤差がｙ軸に、真の数がｘ軸上に示される。

【実施例】

【0065】

〔メソッド例〕
PCR前フラグメントコピー数を推定するための本項目のメソッドは、２つのカテゴリーに入る。項目2.1における第一カテゴリーは、個別観測されたUMIの数をベースにする。これはPCR前後で同じであるため、これらのメソッドはPCRプロセスに依存しない。項目2.2における第二カテゴリーのメソッドは、PCR後のUMIカウント数をベースにする。後者はPCRにより影響を受けるため、これらのメソッドはPCRプロセスについての統計モデルに依存する。項目2.1と2.2のメソッドは、PCR前のUMI分布の情報を必要とする。これは、ヌクレオチドがUMIの各位置に挿入されたことによる仕様からまたはPCR後に入手可能なデータから導き出せる。PCR前のUMI分布はヌクレオチド挿入の頻度により影響を受けるだけでなく、ライゲーションバイアスにも影響されるため、第二のアプローチが好ましいと思われる。該当するメソッドは項目2.3に説明される。
１．準備

【0066】

以下では、UMIの総数をＢにより表し、UMI標識をｂ＝１,.....,Ｂにより表す。UMI ｂのPCR前確率はｐ_bにより表され、そしてPCR前の完全UMI分布は

【数3】

により表される。更に総プール数はＦにより表され、プール標識はｆ＝1,....,Ｆにより表される。プールｆ中のUMI ｂのPCR前カウント数はｋ_bにより表され、そしてプールｆ中の全UMIのPCR前カウント数の完全セットは

【数4】

と記載される。同様に、プールｆ中のUMI ｂのPCR後カウント数はｎ_bにより表され、そして全UMIのPCR後カウント数の完全セットは

【数5】

により表される。

【0067】

更に、

【数6】

を定める。

【0068】

よって、ＫとＮは、PCR前と後のプールｆのサイズである。PCR前とPCR後カウントの当該表記法では異なるプール間を識別しないことに注意されたい。下記では、特に断らない限り、両メソッドが単一プールに適用されることが常に理解されるだろう。

【0069】

２．１ 観測される個別のUMIの数からの推定
2.1.1 観測される個別のUMIのカウント（計数）
UMIカウント数

【数7】

のPCR後セットからPCR前フラグメントコピー数Ｋを推測するために最も単純で最も広く用いられる方法は、

【数8】

で観測される個別のUMIをカウントすることである。この数は、

【数9】

で示されるだろう。次いでＫの推定値は

【数10】

により与えられる。

【0070】

このメソッドは、２つの異なるフラグメントコピーへの同一UMIのライゲーションが全く起こりそうにない場合にのみ、妥当な推定値を生じる。これは、ＢがＫよりはるかに数で勝ること、およびUMIが均等分布（一様分布）することを要求する。

【数11】

がＢに近づくのに十分なほどＫが大きい場合、(1)は真のＫを大幅に過小評価するだろう。
2.1.2 観測される個別のUMIの度数補正カウント

【0071】

観測される個別のUMIの数を単純にカウントすることの代わりに、複製前UMI分布

【数12】

を用いる多峰型（マルチモーダル）分布から導出された個別要素の期待数によってＫを推定することがより妥当と思われる。このメソッドは、頻度補正UMIカウントと称される。ＤＵが個別のUMIの数のランダム変数を示すならば、UMIは

【数13】

となる。

【0072】

上記式は、それがPCR前フラグメントコピーの真の数からのUMIカウント数の偏差を推測するために使われる均等分布UMIについては、参考文献[1]および[2]においても示されている。ここで当該メソッドは

【数14】

であると要求することによってＫを推測する。

【0073】

式（3）は、一般的な

【数15】

については、全く閉鎖解を持たないが、次のように反復して式を解くことができる：

【数16】

【数17】

【0074】

この反復手順は、(3)を解く(5)の不動点に収束する。よって、Ｋ_iは式(3)のＫ_estに収束する。前と同様、このメソッドは、

【数18】

がＢに近づかないためにＮが十分に小さいことを要求する。

【数19】

である極端な事例では、このメソッドは収束しないだろう。

【0075】

２．２ UMIカウントからの推定
前の項目と比較して、次のメソッドはPCR後のUMIカウント数の完全ベクトル

【数20】

を利用する。

【数21】

は複製により影響を受けるので、このメソッドは、PCR用の計算モデルを必要とする。本発明者らは、PCR前に

【数22】

回およびPCR後に

【数23】

回、次の通り、生成するUMIの確率を因数分解する：

【数24】

【0076】

ここで、

【数25】

は、複製効率ｐ_dを有する複製モデルである。後者は、場合により、λとも称される。また

【数26】

は、結果Ｂ、結果の確率

【数27】

および試行Ｋで多峰型（マルチモーダル）分布すると仮定される。すなわち下記となる。

【数28】

【0077】

本発明者らは、本項目で３つの複製モデルを展開する。ポアソンおよび二項分布モデルは一世代の子孫を生成し、一方でガルトン・ワトソン法モデルは複数世代を生成する。後者はPCRプロセスのより正確な記述法であり、従って対応するメソッドは実データ上に良好な結果を生じるはずである。元のフラグメントコピー数Ｋと効率ｐ_dは、理論上、確率を最大化することにより推定することができる：

【数29】

【0078】

しかしながら、(8)中の被加数の数は、Ｋの値が大きい場合には禁止となる。代わりに、本発明では、Ｋとｐ_dは

【数30】

のいずれかを最大化することにより推定される。

【0079】

ここで、

【数31】

は、UMIのセットがｂのみから成る時、PCR後にｎ_b回、UMI ｂおよび別の

【数32】

をそれぞれ確率ｐ_bおよび１－ｐ_bで観測する確率である。（10）における条件付き確率は、更にＮにも依存する。よって、方程式（9）および（10）は、二元UMIセットについてのPCR後UMIカウント数を各々モデル化したＢ非依存性（独立）確率分布の積（総乗）により、

【数33】

を置換する。表記上の簡便さのため、式（9）および（10）中の確率から１－ｐ_bを除き、以降はｑ_b＝１－ｐ_bと設定する。

【0080】

2.2.1 ポアソン複製
ポアソンモデルでは、ｋ個のフラグメントから複製効率λでｎ個のフラグメントを生成する確率は、次のポアソン確率関数により与えられる：

【数34】

【0081】

この場合、（10）の条件付き確率

【数35】

が

【数36】

により与えられる、

【数37】

を有すると証明することができる。

【0082】

（12）と（13）は効率λに依存しないため、このパラメータは以降の表記から除かれるだろう。各UMI bについて、（12）と（13）を正規分布の平均および分散として使用する。すると、Ｋを最大化する（10）は、次の二次多項式の正根として見ることができる。

【数38】

【数39】

ここで

【数40】

【0083】

正規分布によるそれの分散および確率ｐ(ｎ|Ｎ,Ｋ,ｐ_b）のモデル化のため、このメソッドはポアソン－正規またはＰ－正規分布メソッドと呼称される。Ｐ－正規分布メソッドは、複製率λに非依存性であり、（14）と（15）が効率的に計算できるために非常に迅速である。
2.2.2 二項分布複製

【0084】

上記項目でのポアソン複製は、唯１個から無限数のフラグメントを生成できる。しかしながら、実際には、有限数のフラグメントしか生成できない。従って、ポアソンモデルの代わりとなる賢明な置換は二項分布モデルである。後者の場合、複製確率は

【数41】

により与えられる。ここでＭは、１個から生成できるフラグメントの最大数であり、ｐ_dは１フラグメントを複製する１回試行における成功確率である。前と同様、ｑ_d＝１－ｐ_dと定めた。各フラグメントが各PCRサイクルで１回複製されるならば、ｃ回サイクル後は

【数42】

である。

【0085】

更に、１回サイクル後に１フラグメントが複製される確率がｐ_d(1)である場合、ｃ回サイクル後の（16）の成功確率ｐ_d(c)は、

【数43】

により与えられる。Ｍを増加させると確率（16）はλ＝Ｍｐ_dで（11）に収束する。従って、漸近的に、ポアソン複製モデルと二項分布複製モデルは等価である。前と同様、確率

【数44】

の平均は（12）により与えられる。他方、分散

【数45】

は（19）により与えられる：

【数46】

【0086】

式（19）は同様に効率ｐ_dに依存せず、M→∞にすると（13）に収束する。前と同様、Kを最大化する（10）は、（12）と（19）を（10）の正規分布の平均と分散としてそれぞれ使用すると、二次多項式の根として見つけることができる。この場合、多項式は

【数47】

により示され、ここで、Ｖは前と同様に（15）により定義される。Ｋを推定するこのメソッドを、二項－正規またはＢ－正規メソッドと呼ぶことにする。

【0087】

2.2.3 ガルトン・ワトソン複製
PCRの反復性質は、ポアソンまたは二項分布複製モデルよりもむしろ分岐プロセスによってより適切にモデル化される。よって、本発明者らはガルトン・ワトソンプロセスを選択した。本項目で展開するメソッドには、従って、ガルトン・ワトソンまたはGWの接頭辞が付けられる。子孫分布は

【数48】

により示される。本モデルでは、各フラグメントが確率ｐ_dで１サイクルの間に新たな子孫を生みだすことができる。ｋ個のフラグメントが１サイクルの始点に認められる場合、該サイクルの終点にｎ個のフラグメントが検出される確率は二項分布、すなわち

【数49】

である。

【0088】

サイクル数ｃが１より大きい場合、確率

【数50】

は、常母関数とも称されるそれの確率母関数

【数51】

によって最もよく記述される。例えば、参考文献[3]の第10章（Chapter 10）を参照、

【数52】

ここで、

【数53】

はｃ倍組成であり、そして

【数54】

は関数ｆ(ｘ)の第ｋ乗である。これを用いて、次のように書くことができる

【数55】

【0089】

下記において、確率

【数56】

のモデルとして単峰型分布および多峰型分布を調べ、

【数57】

と

【数58】

という事実を利用する。

【0090】

Kを最大化する（9）の検索を単純化するために、更に、ｐ_dが

【数59】

を満たすことを要求する。

【0091】

単峰型（ユニモーダル）分布
（25）と（24）におけるガルトン・ワトソン複製モデルの平均と分散は、それぞれ

【数60】

により示される。

【0092】

本発明者らは、パラメータ分布の平均と分散を式（28）と（29）に設定することにより、ガルトン・ワトソン複製プロセスの単峰型分布モデルを定義する。本発明者らの実験では、平均と分散に関するパラメータ化がそれぞれ

【数61】

により与えられるガルトン・ワトソン複製プロセスの単峰型分布モデルを使用した。

【0093】

多峰型（マルチモーダル）分布
最後の項目のモデルを使用した実験において、ガルトン・ワトソン複製モデルの分布が単一モードでは適切に記述されないことが判明した。従って、本項目では、（25）の混合形態の各成分を別々にモデル化する。PCRをｎ回実施後にｂが検出される確率は、それがｋ回前に観測されると仮定すると、次の平均と分散を有するという事実を利用する。

【数62】

【0094】

デルタ分布と単峰型分布の混合形態
式（33）および（34）は、ｋ＝０の場合、分布ｐ(ｎ|ｋ,ｐ_d,ｃ）はクロネッカーのδ₀に等しいことを示す。よって、ｐ(ｎ|Ｋ,ｐd,ｃ）は、小さなＫの場合ｎ＝０に顕著なピークを有する。ｐ(ｎ|Ｋ,ｐ_d,ｃ）を一様分布でモデル化するのではなく、２つの事例ｎ＝０およびｎ＞０を別々にモデル化するのが妥当である。これはｐ(ｎ|Ｋ,ｐ_b,ｐ_d,ｃ）を混合形態

【数63】

として記述することができる。

【0095】

重みω₁,ω₂、平均μ₁,μ₂、および標準偏差σ₁,σ₂を有する二成分からなる混合分布の場合、該混合分布の全体平均μおよび標準偏差σについての式に従って

【数64】

を得る。よって、μ₂,σ₂はμ,σ,μ₁,σ₁から導出することができる。第一成分がｎ＝０に相当しそして第二成分がｎ≠に相当する場合、ω₁は

【数65】

により与えられる。

【0096】

μ₁およびσ₁は

【数66】

であるので、μ₂については

【数67】

そしてσ₂ ²については

【数68】

と書ける。

【0097】

上記と同様、その平均と分散を（39）と（40）に設定することにより、

【数69】

をパラメトリック確率分布を用いて概算することができる。同様に、この目的に正規分布（30）、ガンマ分布（31）および負の二項分布（32）を使用した。

【0098】

多成分混合分布
例えデルタ分布と単峰型分布の混合形態が小さいＫと大きなＫの状況を正しく記述しているとしても、本発明者らの実験において中間域Ｋについて系統的バイアスを発見した。これは、この領域では総和（25）が小さいｋ項に支配され、そのため負の二項分布またはガンマ分布によるｐ(ｎ|ｋ,ｐ_d,ｃ）の近似値が不正確であるという事実によるものである。これは図６に見ることができ、効率ｐ_dの増加とともに近似の不正確さが増加することも示している。事実、１に近いｐ_dについては、これらの図は、全体の形状がより小さいスケールで繰り返される、フラクタル様構造を示す。（25）の混合成分の複雑構造に適応させるために、ｐ(ｎ|Ｋ,ｐ_b,ｐ_d,ｃ）が（33）と（34）によってその平均と分散が示される３つ以上の分布の混合物である、代替モデルを研究した。ｋ＝０の場合には、再びクロネッカーのδ₀を利用する。ｋ＞０の場合、最初に負の二項分布を適合させた。しかしながら、これは小さなｋにはうまく適合しないことが判明した。従って、ｋを小さい値、中程度の値、大きい値の各範囲に分割することに頼った。小さい値の場合には、ｐ(ｎ|ｋ,ｐ_d,ｃ）は（24）の助けを借りて各ｋについて正確に算出し、中間の値の場合には、各ｋについて（33）と（34）により与えられる平均と分散を有する正規分布を利用し、そして大きいｋの部分和には単一の負の二項分布を利用する。ｋの範囲の区切りには固定境界を選択する。しかしながら、原則的には、それらの境界はモデル分布と真の分布の間の誤差を最小化するように適宜選択することができる。本発明の実験では、小さいｋはｋ＝０からｋ＝15までに及び、中間のｋはｋ＝16からｋ＝49までに及び、そして大きいｋはｋ＝50から始まる範囲内であった。

【0099】

確率分布の切り捨て
ｐ(ｎ|ｋ,ｐd,ｃ）が、ｎ＜ｋのときｆ(ｎ)＞０を有する分布ｆ(ｎ）によりモデル化されるならば、ｆ(ｎ)は範囲ｎ≧ｋに制限されなければならないことに留意すべきである。これは、

【数70】

（ここでｃdｆ_f(ｎ)はｆ(ｎ）の累積関数である）を意味する。これは、もちろん、ｐ(ｎ|ｋ,ｐ_d,ｃ）の平均と分散が無制限のｆ(ｎ）のものと合致しないことも意味する。一般に、分布ｐdｆ（ｎ）が別の分布ｆ(ｎ）の範囲ｎ≧ｋへの制限によって与えられるならば、

【数71】

でｐdｆ(ｎ)＝ｃｆ(ｎ）となり、無制限ｆ(ｎ）の平均と分散は

【数72】

に等しくならなければならない。

【0100】

（42）と（43）の計算は、ｆ(ｎ）が負でない整数（非負整数）上に、例えば負の二項分布またはガンマ分布に定義されれば簡単である。しかしながら、切り捨てはマイナーな役割を果たすことがわかり、大きいｋの場合に１つの負の二項分布の多成分混合物においてのみそれを利用する。

【0101】

2.3 データからのUMI度数の推定
（12）はポアソン、二項分布およびガルトン・ワトソン複製モデルに当てはまるので、これはそれらのモデルから作成されたデータに関して、ｐ_bはｎ_bの分布から推定できることを示す。次のように(12）を書き直したものは

【数73】

【数74】

の期待値がNとKに依存しないことを示す。従って、ｐ_bを推定するための１つの方策は、同一フラグメントの全プールＦに渡り、

【数75】

の平均を計算すること、すなわち

【数76】

であり、これは、ｂに関する確率分布を与える。（44）の期待値はＮおよびＫに対して独立であるが、

【数77】

の分散はそうでない。異常値の存在を避けるため、従って、同様にサイズ指定されたＮの十分に大きいセットに渡って

【数78】

の平均を計算することが意味をなすに違いない。別の代替法は、

【数79】

によりｐ_bを推定するものである。これは、ｐ_bの推定値に対する小さなプールの影響を減少させる。

【0102】

３．実験
3.1 合成データの作成
当該実験では、均等および不均等分散しているＢ＝64およびＢ＝256 UMIを使用した。不均等データについては、UMIの各位置でのヌクレオチド分布は、ディリクレ分布から、すなわち

【数80】

から抽出した。ここでｐ_A、ｐ_c、ｐ_G、ｐ_Tは、ヌクレオチドＡ,Ｃ,ＧおよびＴを観測する確率である。全体的に見れば、UMI bの確率は

【数81】

により与えられ、ここでnt(b,i)はバーコードｂのｉ番目の位置のヌクレオチドであり、その積（総乗）はｂの全位置に渡る。式（47）においてα＝５と設定した。Ｂ＝64とＢ＝256についての対応する確率密度関数ｐ_bは、図７に与えられる。この図面は、１／Ｂの周辺のｐ_bのばらつき（垂線）とそれらの確率密度関数の非対称形を示す。よって、抽出されたUMI分布

【数82】

は均等分布からある程度かけ離れている。本発明者らは、10から10000までのＫについて本発明方法を試験した。各Ｋについて、50回から100回の間で

【数83】

を抽出し、 そして各標本

【数84】

について、多項分布（7）に従って

【数85】

を得た。

【0103】

次に、ポアソン複製では（11）により示され、二項分布複製では（16）により示され、そしてガルトン・ワトソン分布複製では（24）により示される、（6）の確率

【数86】

から

【数87】

を抽出した。ガルトン・ワトソンモデルでは効率

【数88】

を0.5、0.7、0.8および0.9に設定し、サイクル数を15に設定した。ポアソンモデルと二項分布モデルには、λと

【数89】

は式

【数90】

（ここでＭ＝２^c－１である）
により、ｐ_{d_GW}から導出される。以下では、λとｐ_{d_Bin}よりもｐ_{d_GW}を参照し、ｐ_{d_GW}を単純にｐ_dと記すことにする。

【0104】

3.2 結果
Ｂ＝64とＢ＝256 UMIについてのUMIカウントの結果は図８に見ることができるが、それらの頻度補正したUMIカウントは図９に与えられる。別のメソッドに比較して、UMIカウントの両バージョンは、複製モデルに依存せず（独立であり）、従って複製効率に影響されない。この理由により、効率0.5を有する二項分布モデルの結果のみを含めた。他の複製モデルと効率の結果は同等である。図８は、UMIカウントが、UMIの数ＢがＫよりも相当大きい場合にのみ、良好なＫの推定値を生じることを示す。Ｂ＝64に関しては、推定値がＫ＝10についてのみ妥当な正確さを有する。Ｂ＝256に関しては、Ｋ＝100の推定も妥当である。しかしながら、より大きいＫの値には、UMIカウントはＫを少なくとも50％過小評価する。図８は更に、UMIカウントが不均等データに関してわずかに低い精度を生じる。その一方、図９に記載の頻度補正したカウントは、より広範な領域のＫについて偏りのない（バイアスフリーの）より正確な結果をもたらす。しかしながら、大きなＫの場合には、観測される個別のUMIの数

【数91】

がＢに近づき、このメソッドはＫを収束させることに失敗するかまたはＫを大幅に過小評価してしまう。不均等なバーコード分布は、当メソッドが収束する領域を拡張する。これは、不均等バーコード分布では全てのUMIを観測することがより困難であり、よって

【数92】

がより高位のＫについてのみＢに近づくという事実によるものである。他方で、不均等UMI分布は、小さいＫ値について精度の低下を引き起こす。

【0105】

（14）および（15）の二次方程式の正根として定義される、ポアソン－正規分布メソッドの結果は、Ｂ＝64とＢ＝256についてそれぞれ図１０と１１に与えられる。それらの結果は、Ｐ－正規メソッドの精度が効率ｐ_dにほとんど影響されない（独立性である）ことを示す。しかしながら、不均等性は、小さいＫに対して負の効果を有するようである。頻度補正したカウントに比較して、Ｐ－正規分布メソッドは、全領域のＫに渡って良好な結果を与えた。図１２および１３は、二項－正規分布メソッドの結果を示す。これは、ポアソン－正規分布メソッドと二項－正規分布メソッドの間には非常にわずかな相違しかないことを示す。項目2.2.2に記載の通り、これはポアソン分布モデルと二項分布モデルの残金的等価性のためである。ポアソン－正規分布メソッドと二項－正規分布メソッドとの間のわずかな相違は、小さいＫについてのみ観察され得る。

【0106】

図１４と１５は、全ＧＷメソッドの中で最高の精度を有したＧＷ－多成分混合メソッドの結果を含む。図１４と１５は、効率ｐ_dが増加するとともに精度がわずかに増加することを示す。他方で、不均等性（非一様性）は結果にほとんど影響を与えない。ＧＷ－多成分混合メソッドは、全領域のＫに対して有益な結果をもたらし、特に小さなＫについてポアソン－正規および二項－正規分布メソッドよりも正確である。ＧＷ－多成分混合メソッドのＢ＝256 UMIおよびかなり高い精度は、他の全てのメソッドを大幅に上回る。このメソッドの優越性を考慮すると、ポアソンおよび二項分布複製モデルについても多成分混合を用いて

【数93】

をモデル化することは価値があるかもしれない。

【0107】

要約すれば、本発明者らの実験は、当該メソッドが複製前フラグメントの数Ｋについての正確な推定値を与えることを証明する。これらの推定値は、UMIの複製効率ｐ_dおよび複製前分布

【数94】

に大部分独立である。頻度補正したUMIカウントとは異なり、本発明メソッドは全て、調査したPCR前フラグメント番号Ｋの全領域について結果を与える。

【図1】

【図2】

【図3】

【図4】

【図5】

【図6】

【図7】

【図8】

【図9】

【図10】

【図11】

【図12】

【図13】

【図14】

【図15】