特許 7423101 細胞集団の処理方法および細胞集団に含まれる遺伝子の分析方法

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

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

(11)【特許番号】

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

(45)【発行日】2024-01-29

(54)【発明の名称】細胞集団の処理方法および細胞集団に含まれる遺伝子の分析方法

(51)【国際特許分類】

   C12Q 1/6869 20180101AFI20240122BHJP

   C12Q 1/686 20180101ALI20240122BHJP

   C12N 15/09 20060101ALN20240122BHJP

   C12Q 1/6806 20180101ALN20240122BHJP

【ＦＩ】

C12Q1/6869 Z ZNA

C12Q1/686 Z

C12N15/09 Z

C12Q1/6806 Z

【請求項の数】 16

(21)【出願番号】P 2022526658

(86)(22)【出願日】2021-05-28

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

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

(87)【国際公開日】2021-12-02

【審査請求日】2022-10-19

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

(32)【優先日】2020-05-29

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

【国等の委託研究の成果に係る記載事項】（出願人による申告）平成３０年度、生体システム理解・医科学応用を実現する１細胞核酸計測技術の開発、産業技術力強化法第１７条の適用を受ける特許出願

【早期審査対象出願】

(73)【特許権者】

【識別番号】503359821

【氏名又は名称】国立研究開発法人理化学研究所

(74)【代理人】

【識別番号】230104019

【弁護士】

【氏名又は名称】大野 聖二

(74)【代理人】

【識別番号】100173185

【弁理士】

【氏名又は名称】森田 裕

(72)【発明者】

【氏名】城口 克之

(72)【発明者】

【氏名】ジン ジャンシ

(72)【発明者】

【氏名】山本 れいこ

【審査官】上村 直子

(56)【参考文献】

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

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

【文献】米国特許出願公開第２０１９／００２５２９９（ＵＳ，Ａ１）

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

【文献】米国特許出願公開第２０１９／０００２９７４（ＵＳ，Ａ１）

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

【文献】米国特許出願公開第２０１９／０３００９６８（ＵＳ，Ａ１）

【文献】Nature Communications，2020年，11:866，pp.1-9，Published online: 13 February 2020

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

Ｃ１２Ｑ １／００－１／７０

Ｃ１２Ｍ １／００－１／４２

ＪＳＴＰｌｕｓ／ＪＭＥＤＰｌｕｓ／ＪＳＴ７５８０（ＪＤｒｅａｍＩＩＩ）

ＣＡｐｌｕｓ／ＭＥＤＬＩＮＥ／ＥＭＢＡＳＥ／ＢＩＯＳＩＳ（ＳＴＮ）

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

【請求項1】

細胞集団に含まれる所定の遺伝子の塩基配列を分析するための方法であって、
（Ａ）単離された細胞集団と固有のバーコード配列をそれぞれ有する細胞バーコードと核酸増幅用のプライマーおよび遺伝子増幅用試薬を含む細胞の分散液から、液滴集団であって、水性の液滴を含み、前記液滴の少なくとも一部のそれぞれが、１つの細胞と１分子の細胞バーコードとを含む液滴集団を得ることと、ここで、前記液滴集団は、細胞バーコードの増幅産物と前記細胞集団の細胞内の所定の遺伝子それぞれの増幅産物を得るために必要な核酸増幅用のプライマーおよび遺伝子増幅用試薬を含み、
（Ｂ）得られたそれぞれの液滴中で、細胞バーコードの増幅産物と所定の遺伝子それぞれの増幅産物を得て、さらに、細胞バーコードと所定の遺伝子の全部または一部の塩基配列を含む連結物を得ることと、得られた連結物を液滴から水溶液中に回収して、得られた連結物をシークエンスして所定の遺伝子の塩基配列と細胞バーコードの塩基配列を決定することと、
（Ｃ－１）決定された細胞バーコードの塩基配列に基づいて、決定された塩基配列をクラスタリングして、複数の第一のクラスターを得ることと、
を含む、
方法。

【請求項2】

前記（Ｂ）において、細胞バーコードの増幅産物は、第一のプライマーに由来する第一の領域を有し、所定の遺伝子の増幅産物は、第二プライマーに由来する第二の領域を有し、第一の領域と第二の領域は、互いにハイブリダイズ可能な相補的な配列部分を有し、前記第一のプライマーおよび第二のプライマーはそれぞれ、１以上のタグ分子を連結しており、当該タグ分子は、前記連結物には含まれず、かつ、
前記（Ｂ）において、水溶液中に回収された連結物から、タグ分子を有する増幅産物を当該タグ分子に親和性を有する分子を担持したカラムまたはビーズを用いて除去することをさらに含む、請求項１に記載の方法。

【請求項3】

（Ｄ－１）得られた第一のクラスターの数から細胞集団に含まれる細胞の数または特定の所定の遺伝子を有する細胞の数を推定すること
をさらに含む、請求項１または２に記載の方法。

【請求項4】

（Ｃ－２）決定された所定の遺伝子の塩基配列に基づいて、決定された塩基配列をクラスタリングして、複数の第二のクラスターを得ること
をさらに含む、請求項１～３のいずれか一項に記載の方法。

【請求項5】

（Ｄ－２）得られた第二のクラスターの数から細胞集団に含まれる細胞の種類の数を推定することをさらに含む、請求項４に記載の方法。

【請求項6】

（Ｄ－３）得られた細胞バーコードの塩基配列と所定の遺伝子の塩基配列の組合せの情報に基づいて、少なくとも１つのある第二のクラスターに分類された所定の遺伝子の塩基配列と連結している細胞バーコードの塩基配列から当該所定の遺伝子の塩基配列が分類された第一のクラスターを決定し、当該細胞バーコードが分類された第一のクラスターの数から、当該第二のクラスターに分類された細胞の数を推定すること
をさらに含む、請求項４に記載の方法。

【請求項7】

（Ｃ－４）同一の第一のクラスターに分類された配列が異なる第二のクラスターに分類される場合、当該第二のクラスターを同一の細胞ベースの操作上分類単位（ｃＯＴＵ）に分類することをさらに含む、請求項４に記載の方法。

【請求項8】

（Ｅ）第一の細胞集団と、第一の細胞集団とは異なる第二の細胞集団のそれぞれに関して、細胞集団に含まれる（ｉ）ｃＯＴＵの数および／または（ｉｉ）特定のｃＯＴＵに含まれる細胞の数を推定し、第一の細胞集団に関して推定された（ｉ）ｃＯＴＵの数および／または（ｉｉ）特定のｃＯＴＵに含まれる細胞の数を、第二の細胞集団に関して推定された（ｉ）ｃＯＴＵの数および／または（ｉｉ）特定のｃＯＴＵに含まれる細胞の数とそれぞれ比較することをさらに含む、請求項７に記載の方法。

【請求項9】

（Ｆ）第一の細胞集団に関して推定された（ｉ）ｃＯＴＵの数および（ｉｉ'）特定のｃＯＴＵに含まれる細胞の数と、第二の細胞集団に関して推定された（ｉ）ｃＯＴＵの数および（ｉｉ'）特定のｃＯＴＵに含まれる細胞の数を比較することを含む、請求項８に記載の方法。

【請求項10】

細胞集団が、微生物叢である、請求項１～９のいずれか一項に記載の方法。

【請求項11】

微生物叢が、体内または体表の微生物叢である、請求項１０に記載の方法。

【請求項12】

微生物叢が、消化管内の微生物叢である、請求項１０に記載の方法。

【請求項13】

第一の細胞集団と第二の細胞集団が、同一対象の異なる部位から取得された微生物叢である、請求項８または９に記載の方法。

【請求項14】

第一の細胞集団と第二の細胞集団が、異なる対象の同一の部位から取得された微生物叢である、請求項８または９に記載の方法。

【請求項15】

第一の細胞集団と第二の細胞集団が、同一対象の同一の部位から異なる時間に取得された微生物叢である、請求項８または９に記載の方法。

【請求項16】

細胞集団が、未知の細胞を含む、請求項１～１５のいずれか一項に記載の方法。

【発明の詳細な説明】

【技術分野】

【0001】

本発明は、細胞集団の処理方法および細胞集団に含まれる遺伝子の分析方法に関する。

【背景技術】

【0002】

共生微生物叢の構成がどのように宿主の健康^１，２に寄与するかを本質的に理解するために、細胞が微生物叢^３－５の基本的で物理的な単位であるから、単純に微生物叢を細胞レベルで定義すべきである。しかし、現在の最新の技術では困難である^６－８。

【0003】

微生物相と宿主の相互作用は宿主のホメオスタシスおよび多くの疾患^{１３－１６}と関連している。微生物相と宿主の相互作用のメカニズムをさらにかつ統合的に理解するために、微生物相を研究するだけでなく、微生物相と宿主の両方についてメタボロミクスおよび／またはトランスクリプトミクスのような他の分析と微生物相の組成分析を結びつけることが重要である^５。この目的のためには、一般に使用可能な単位、例えば、重量あたりの細胞数および／または体積あたりの分子数に基づく濃度測定が必要である。この点に関して、細胞内に存在する核酸の分子数をカウントする技術が開発されている（特許文献１～３）。このカウント技術では、分子一つ一つに固有の核酸配列（バーコード）を付与して、バーコードの種類数をカウントすることによって分子数を推定する。また、特許文献１～３において核酸の増幅中のエラーやシークエンス時の読み取りエラーにより、分子のカウント数に誤差が生じ得る。この誤差を低減する技術も開発されている（特許文献４）。特許文献４では、核酸の増幅中のエラーやシークエンス時の読み取りエラーの性質を考慮した、エラー除去とカウント数の是正に関する方法が提案されている。しかし、現行の技術^６－８では、細胞レベルでの微生物相組成の測定は困難であった。さらに、微生物叢は多数の細菌種による膨大な数の細菌からなる^１７。しかしながら、高い分類学的分解能を有するハイスループットの細胞定量法はこれまでに開発されていない。

【0004】

次世代シークエンシング技術を用いた１６Ｓ ｒＲＮＡ遺伝子アンプリコンシークエンシングに基づいたハイスループット法が、細菌の多様性の研究に貢献した^{２２，２３}。しかし、従来の方法は、精製バルク細菌ゲノムから１６Ｓ ｒＲＮＡ遺伝子を増幅し、増幅した分子の数を測定するので、基本的には以下の限界がある。１）異なる種はゲノム上に１６Ｓ ｒＲＮＡ遺伝子の異なるコピー数を有し、大部分の種のコピー数は不明であるため、細胞数を測定し、異なる種の細胞数を比較するのは困難である。２）１６Ｓ ｒＲＮＡ配列の同定は、配列決定および増幅エラーのために正確ではなく、低分類学的分解能をもたらす。実際、配列決定エラーは、分子バーコード^{２４－２６}を用いて修正されたが、主にキメラ生成に由来する増幅エラーを十分に除去できない^２７。

【先行技術文献】

【特許文献】

【0005】

【文献】ＵＳ９２６０７５３Ｂ

【文献】ＵＳ１０２８７６３０Ｂ

【文献】ＵＳ１０５８４３８２Ｂ

【文献】ＷＯ２０１８／２３５９３８

【発明の概要】

【0006】

本発明は、細胞集団の処理方法および細胞集団に含まれる遺伝子の分析方法を提供する。

【0007】

本発明者らは、細菌微生物叢の細胞型と各細胞型に対する細胞濃度をハイスループット法で定量する新しい方法を開発した。本発明者らはまた、１つの細胞に解析対象の遺伝子が重複して存在することに対応した、細胞ベースの操作上分類単位（ｃＯＴＵ）に基づいて、解析対象の遺伝子群を分類することによって、遺伝子重複を有する未知の細胞（例えば、微生物）の細分類を可能とし、その数を推定する方法を見出した。

【0008】

本発明によれば以下の発明が提供される。
［１］細胞集団の処理方法であって、
（Ａ）単離された細胞集団を含む細胞の分散液から、液滴集団であって、水性の液滴を含み、前記液滴の少なくとも一部のそれぞれが、１つの細胞と１分子の細胞バーコードとを含む液滴集団を得ること
を含む、方法。
［２］細胞集団に含まれる遺伝子の塩基配列を分析する方法であって、
（Ａ）単離された細胞集団を含む細胞の分散液から、液滴集団であって、水性の液滴を含み、前記液滴の少なくとも一部のそれぞれが、１つの細胞と１分子の細胞バーコードとを含む液滴集団を得ることと、
（Ｂ）得られたそれぞれの液滴中で、細胞バーコードの増幅産物と所定の遺伝子の増幅産物を得て、さらに、細胞バーコードと所定の遺伝子の全部または一部の塩基配列を含む連結物を得ることと、得られた連結物を液滴から水溶液中に回収して、得られた連結物をシークエンスして所定の遺伝子の塩基配列と細胞バーコードの塩基配列を決定することとを含む、方法。
［３］前記（Ｂ）において、細胞バーコードの増幅産物は、第一のプライマーに由来する第一の領域を有し、所定の遺伝子の増幅産物は、第二プライマーに由来する第二の領域を有し、第一の領域と第二の領域は、互いにハイブリダイズ可能な相補的な配列部分を有し、前記第一のプライマーおよび第二のプライマーはそれぞれ、１以上のタグ分子を連結しており、当該タグ分子は、前記連結物には含まれず、かつ、
前記（Ｂ）において、水溶液中に回収された連結物から、タグ分子を有する増幅産物を当該タグ分子に親和性を有する分子を担持したカラムまたはビーズを用いて除去することをさらに含む、上記［２］に記載の方法。
［４］（Ｃ－１）決定された細胞バーコードの塩基配列に基づいて、決定された塩基配列をクラスタリングして、複数の第一のクラスターを得ること
をさらに含む、上記［２］または［３］に記載の方法。
［５］（Ｄ－１）得られた第一のクラスターの数から細胞集団に含まれる細胞の数または特定の所定の遺伝子を有する細胞の数を推定すること
をさらに含む、上記［４］に記載の方法。
［６］（Ｃ－２）決定された所定の遺伝子の塩基配列に基づいて、決定された塩基配列をクラスタリングして、複数の第二のクラスターを得ること
をさらに含む、上記［２］または［３］に記載の方法。
［７］（Ｄ－２）得られた第二のクラスターの数から細胞集団に含まれる細胞の種類の数を推定することをさらに含む、上記［６］に記載の方法。
［８］（Ｃ－３）決定された細胞バーコードの塩基配列に基づいて、決定された塩基配列をクラスタリングして、複数の第一のクラスターを得ることと、決定された所定の遺伝子の塩基配列に基づいて、決定された塩基配列をクラスタリングして、複数の第二のクラスターを得ること
をさらに含む、上記［２］または［３］に記載の方法。
［９］（Ｄ－３）得られた細胞バーコードの塩基配列と所定の遺伝子の塩基配列の組合せの情報に基づいて、少なくとも１つのある第二のクラスターに分類された所定の遺伝子の塩基配列と連結している細胞バーコードの塩基配列から当該所定の遺伝子の塩基配列が分類された第一のクラスターを決定し、当該細胞バーコードが分類された第一のクラスターの数から、当該第二のクラスターに分類された細胞の数を推定すること
をさらに含む、上記［８］に記載の方法。
［１０］（Ｃ－４）同一の第一のクラスターに分類された配列が異なる第二のクラスターに分類される場合、当該第二のクラスターを同一の細胞ベースの操作上分類単位（ｃＯＴＵ）に分類することをさらに含む、上記［８］に記載の方法。
［１１］（Ｅ）第一の細胞集団と、第一の細胞集団とは異なる第二の細胞集団のそれぞれに関して、細胞集団に含まれる（ｉ）ｃＯＴＵの数および／または（ｉｉ）特定のｃＯＴＵに含まれる細胞の数を推定し、第一の細胞集団に関して推定された（ｉ）ｃＯＴＵの数および／または（ｉｉ）特定のｃＯＴＵに含まれる細胞の数を、第二の細胞集団に関して推定された（ｉ）ｃＯＴＵの数および／または（ｉｉ）特定のｃＯＴＵに含まれる細胞の数とそれぞれ比較することをさらに含む、上記［１０］に記載の方法。
［１２］（Ｆ）第一の細胞集団に関して推定された（ｉ）ｃＯＴＵの数および（ｉｉ’）特定のｃＯＴＵに含まれる細胞の数と、第二の細胞集団に関して推定された（ｉ）ｃＯＴＵの数および（ｉｉ’）特定のｃＯＴＵに含まれる細胞の数を比較することを含む、上記［１１］に記載の方法。
［１３］細胞集団が、微生物叢である、上記［１］～［１２］のいずれかに記載の方法。［１４］微生物叢が、体内または体表の微生物叢である、上記［１３］に記載の方法。
［１５］微生物叢が、消化管内の微生物叢である、上記［１３］に記載の方法。
［１６］第一の細胞集団と第二の細胞集団が、同一対象の異なる部位から取得された微生物叢である、上記［１１］または［１２］に記載の方法。
［１７］第一の細胞集団と第二の細胞集団が、異なる対象の同一の部位から取得された微生物叢である、上記［１１］または［１２］に記載の方法。
［１８］第一の細胞集団と第二の細胞集団が、同一対象の同一の部位から異なる時間に取得された微生物叢である、上記［１１］または［１２］に記載の方法。
［１９］細胞集団が、未知の細胞を含む、上記［１］～［１８］のいずれかに記載の方法。

【図面の簡単な説明】

【0009】

【図1】ＢａｒＢＩＱとその有効性。ａ、ＢａｒＢＩＱの概略図。試料を溶液に懸濁させた後、細菌の塊を破壊するためにボルテックス処理を行った。細胞バーコード、細胞に固有の塩基（１分子毎に異なる塩基配列）を含み、増幅のためのプライム部位を含むＤＮＡ分子；プライマー、１６Ｓ ｒＲＮＡ遺伝子および細胞バーコードのそれぞれの増幅のための、両方の増幅産物を連結するための、およびシークエンスアダプターの付着のためのＤＮＡプライマー；試薬、ＤＮＡ増幅のための試薬。ライブラリーの生成、精製、およびシークエンスの概略図の詳細は、図５、データ処理の詳細は、図６参照のこと。ｂ、模擬細胞集団によるＢａｒＢＩＱの有効性。編集距離、置換、挿入および欠失の最小数として定義されるレーベンシュタイン距離^２９；Ｓａｎ配列、サンガーシークエンスにより同定された１６Ｓ ｒＲＮＡ配列；ＡＴＣＣ／ＪＣＭ／ＤＳＭ－＜数＞、株ＩＤ；Ａ、Ｂ、またはＣ、各株に対するＳａｎ配列；Ｂａｒ配列－ＭＫ－ＸＸ（０１－１６）、ＢａｒＢＩＱで同定された配列（Ｂａｒ配列）；ＣＯＴＵ－ＭＫ－ＸＸ（０１－１０）、細胞ベースの操作上分類単位（ｃＯＴＵ）；赤いアスタリスク記号、１塩基相違を有するＢａｒ配列；ＯＴＵ－ＲｅｐＳｅｑ－ＭＫ－ＸＸ（０１－１２）、ＯＴＵを表す配列；ｃ、ＢａｒＢＩＱ［Ｃ］_{ＢａｒＢＩＱ}および顕微鏡画像［Ｃ］_顕微鏡（表１および２のデータ）により測定した模擬細胞集団における１０株の細胞濃度の比較。青線：対数目盛りで傾き１が固定されたフィッティングライン；ｒ、Ｐｅａｒｓｏｎ係数、Ｒ２：決定係数。エラーバー、標準偏差（［Ｃ］_{ＢａｒＢＩＱ}についてはｎ＝３、［Ｃ］_顕微鏡についてはｎ＝５）。

【図2】マウス盲腸微生物叢の包括的分析。ａ、マウス盲腸における遠位（ｄｉｓｔ）および近位（ｐｒｏｘ）のサンプリング位置。ｂ、Ｂａｒ配列の配列同一性プロフィール；同一性、各Ｂａｒ配列と、３つの一般的な公共データベースＧｒｅｅｎＧｅｎｅ（ＧＧ）、リボソームデータベースプロジェクト（ＲＤＰ）、およびＳｉｌｖａにおけるその最も近い１６Ｓ ｒＲＮＡ配列間との同一性。Ｔｈｒｅｅ、３つすべてのデータベースの組み合わせ。ｃ、技術的反復間のｃＯＴＵの細胞濃度の比較（他の反復は、図１６）。マゼンタ線、ポアソン分布および総濃度による正規化に基づくサンプリングノイズの理論的信頼区間（９９．９％）；水色線、２倍の変化；青い点、異なる濃度を示すｃＯＴＵ；挿入された数字、青色および灰色のドットの数；Ｍａ、Ｍｂ、Ｍｃ、マウス；ｄｉｓｔおよびｐｒｏｘ、位置；１、２、および３、技術的反復。ｄとｅ、ｃと同じであるが、異なる試料間の比較；異なるｃＯＴＵ（青い点）の最小（ｄ）と最大（ｅ）の数の例。他のサンプルの比較については、図１６。ｆ、各サンプル間の微生物のＢｒａｙ－Ｃｕｒｔｉｓの非類似性。ラベルはｃ～ｅと同じである。

【図3】マウス間における個々のｃＯＴＵの細胞濃度の変動。ａ、マウス３匹（Ｍａ、ＭｂおよびＭｃ）における遠位（赤色実線）および近位（水色破線）のｃＯＴＵ細胞濃度の例。ＣＶ、変動係数。ｂ、Ｃｌｏｓｔｒｉｄｉｕｍ ＸＩＶａ属のｃＯＴＵのＣＶ（検出されたすべての属を図９ａおよび９ｂに示す）。ＣＯＴＵ－ＣＭ－＜数＞、ｃＯＴＵのＩＤ；遠位および近位、位置；エラーバー、サンプリングノイズおよび総濃度測定の技術的誤差を想定したシミュレーションによって得られた各ｃＯＴＵに対するＣＶの９５％信頼区間。

【図4】相関細菌ネットワーク。ａ、ｃＯＴＵペアの豊富度に基づく相関の例。ドット、６サンプル（Ｍａ^ｄｉｓｔ１、Ｍａ^ｐｒｏｘ１、Ｍｂ^ｄｉｓｔ、Ｍｂ^ｐｒｏｘ、Ｍｃ^ｄｉｓｔ、Ｍｃ^ｐｒｏｘ）の細胞濃度（細胞／ｍｇ）；ｒ、Ｐｅａｒｓｏｎ係数。ｂ、強く相関する細菌群（ＳＣＢＧ）の定義。樹状図、定義された距離、１－最小（│ｒ’│）［ｒ’∈（ｒ－ＯＣＩ，ｒ＋ＯＣＩ）］に基づき６つのサンプル全てで一般的に検出される２９６ｃＯＴＵの階層的クラスタリング；赤色の破線、閾値０．６；ヒートマップ、すべてのｃＯＴＵのｒ；ヒートマップの白いギャップ、垂直および水平の両方で閾値０．６より低い枝の分離を示す間隔；下の数字、ＳＣＢＧのＩＤ。全ＳＣＢＧのｃＯＴＵ名とＩＤを有する樹状図を図１７に示す。ｃ、それぞれＦｏｒｃｅ－ｄｉｒｅｃｔｅｄ ｌａｙｏｕｔ^３９で可視化したＳＣＢＧ７とＳＣＢＧ２６のｃＯＴＵネットワーク。ノード、ｃＯＴＵ；ノードサイズ、ａのように６つのサンプル中のｃＯＴＵの平均細胞濃度；エッジ色、端で連結されたｃＯＴＵ間のｒ。全ＳＣＢＧの可視化ネットワークを図１２ａ～ｆに示す。ｄ、Ｆｏｒｃｅ－ｄｉｒｅｃｔｅｄ ｌａｙｏｕｔによって可視化されたＳＣＢＧのネットワーク。エッジ色、２つのＳＣＢＧ間のＳＣＢＧ間相互関係Ｒ_{ｉｎｔｅｒ}。

【図5】ＢａｒＢＩＱにおけるライブラリー生成の配列情報、精製およびシークエンシングを含む概略図。Ｉ、ＩＩ、ＩＩＩ、ＩＶはＰ５－ｉｎｄｅｘ－Ｒ１Ｐ－ｂａｒｃｏｄｅ－Ｒ、Ｂｉｏｔｉｎ－Ｌｉｎｋ－ｂａｒｃｏｄｅ－Ｆ、Ｂｉｏｔｉｎ－ｌｉｎｋ－８０５Ｒ、Ｐ７－Ｒ２Ｐ－３４１Ｆと命名され、設計されたプライマー；Ｉｎｄｅｘ（ＸＸＸＸＸＸＸＸ）は設計された８塩基；Ｂａｒｃｏｄｅ、ランダムおよび固定塩基（他の３種類のバーコードは表３にある）；配列中のＮ、Ａ、Ｃ、Ｇ、またはＴ；Ｉ２、Ｒ１およびＲ２、ＭｉＳｅｑのためのイルミナシークエンシングプライマー；Ｉ１、カスタマイズされたシークエンシングプライマー。

【図6】ＢａｒＢＩＱデータ処理の概略図。黒矢印、処理段階；赤矢印、次の段階のオペランドの記述；バーコード、細胞バーコード；Ｒ１、Ｒ１のリード；Ｉ１、およびＲ２、低品質の末端とプライマー部分がトリミングされたＩ１、およびＲ２のリード；ＢＣｌｕｓｔｅｒ、バーコードによってクラスタリングされたクラスター；ＳＣｌｕｓｔｅｒ、各ＢＣｌｕｓｔｅｒにおける１６Ｓ ｒＲＮＡの配列によってクラスタリングされたサブクラスター；シフトしたＲｅｐＳｅｑ、リードのプライマー部分における挿入または欠失のために生じたＲｅｐＳｅｑ；１つの挿入または欠失のＲｅｐＳｅｑ、トリミング後にリードの残りの部分における１塩基の挿入または欠失のエラーにより生じたＲｅｐＳｅｑ；キメラＲｅｐＳｅｑ、ＰＣＲキメラによって得られ得るＲｅｐＳｅｑ；稀なエラーＲｅｐＳｅｑ、トリミング後にリードの残りの部分における１つのインデル（挿入または欠失）および１つの置換、１つのインデルおよび２つの置換、または２つのインデルのエラーによって生じるＲｅｐＳｅｑ；ＲｅｐＳｅｑ型、ＲｅｐＳｅｑの配列型；低カウントＲｅｐＳｅｑ、少ない数のＢＣｌｕｓｔｅｒにおいて検出されるＲｅｐＳｅｑの型；１塩基エラーＲｅｐＳｅｑ、他のＲｅｐＳｅｑと１塩基の相違を有するＲｅｐＳｅｑの型、前者および後者のＲｅｐＳｅｑの型の間で検出されたＲｅｐＳｅｑの数が閾値よりも小さい；Ｂａｒ配列、ＢａｒＢＩＱで同定された配列；ｃＯＴＵ、細胞ベースの操作上分類単位。

【図7】ＢａｒＢＩＱ測定において、各ｃＯＴＵの絶対細胞濃度および各ｃＯＴＵのサンプリングノイズを計算するために使用される総濃度。ａ、液滴デジタルＰＣＲ（実施例のＢａｒＢＩＱ法の項目を参照）で測定した各試料の総細菌濃度。Ｍａ、Ｍｂ、Ｍｃ、Ｍｄ（Ｍｄはシークエンスされていない）、マウス；ｄｉｓｔとｐｒｏｘ、位置（図２ａ参照）；１、２、および３、技術的反復。エラーバー、標準偏差（ｎ＝５）。ｂ、カウントの平均の関数として、各ｃＯＴＵに対するＭａ^ｄｉｓｔの技術的反復３回におけるカウントのＣＶ^２（ＣＶ、変動係数）；シミュレーション１および２、ならびに理論値をポアソン分布に基づいて得た。ｃ、ｌｏｇ_１０（ＣＶ^２）－ｌｏｇ_１０（ＣＶ_{Ｐｏｉｓｓｏｎ} ^２）の分布；ＣＶ各、ｃＯＴＵのＣＶ；ＣＶ_{Ｐｏｉｓｓｏｎ}、Ｐｏｉｓｓｏｎ分布に基づく理論的ＣＶ。ｄ、Ｍａ^ｄｉｓｔの測定とシミュレーション１の間、およびシミュレーション１とシミュレーション２の間のｌｏｇ_１０（ＣＶ^２）－ｌｏｇ_１０（ＣＶ_{Ｐｏｉｓｓｏｎ} ^２）の分布のＱ－Ｑプロット^４５。ｌｏｇ_１０（ＣＶ^２）－ｌｏｇ_１０（ＣＶ_{Ｐｏｉｓｓｏｎ}２）の分布は測定とシミュレーションの間で同等であり、検出された各ｃＯＴＵについてのノイズは主にサンプリングによることを示唆した。

【図8】マウスＭａにおける各ｃＯＴＵの位置依存性細胞濃度比較。ａ、マウスＭａにおける遠位位置（Ｍａ^ｄｉｓｔ）と近位位置（Ｍａ^ｐｒｏｘ）間の３回の技術的反復に対する各ｃＯＴＵの平均細胞濃度比較；エラーバー、標準偏差（ｎ＝３）；赤い点、ＦＤＲ＜０．０５および平均＞２の倍率変化を示す（ｂ参照）；破線、倍率変化＝２。ｂ、Ｖｏｌｃａｎｏ ｐｌｏｔは、Ｍａにおける遠位および近位位置の間のｃＯＴＵの細胞濃度の相違を示す。ＦＤＲ（偽発見率）は、関数ｔ．ｔｅｓｔ（Ｒパッケージｓｔａｔ）による両側２群ｔ検定（ｎ＝３）を用いて計算された全２４０のｃＯＴＵのｐ値に基づいてＢＨ法を用いる関数ｐ．ａｄｊｕｓｔ（Ｒパッケージｓｔａｔ）により決定された^４６；Ｍａ^ｄｉｓｔ／Ｍａ^ｐｒｏｘ、Ｍａ^ｐｒｏｘにおける平均細胞濃度に対するＭａ^ｄｉｓｔの平均細胞濃度の比；破線、Ｍａ^ｄｉｓｔとＭａ^ｐｒｏｘの総濃度の比。

【図9a】すべてのｃＯＴＵの分類学に対するＣＶ（変動係数）。左、門から属への分類。右、遠位および近位位置における各ｃＯＴＵのＣＶ。ＣＯＴＵ－ＣＭ－＜数＞、ｃＯＴＵＩＤ；誤差バー、総濃度測定のサンプリングノイズと技術的誤差を想定したシミュレーションで得られた９５％信頼区間。

【図9b】同上。

【図10】ネットワーク全体における各ｃＯＴＵの他のｃＯＴＵとの相関特性。上段、与えられたｃＯＴＵとそれ以外のすべてのｃＯＴＵの間の│ｒ│の分布を示し、│ｒ│は絶対ピアソン相関係数である；ｃＯＴＵは、横軸に沿った各ｃＯＴＵの│ｒ│の平均（シアンの線）によって並べられた。下段、相対周波数によって示された各ｃＯＴＵに対する│ｒ│の分布；上図の各行について、数値は、それらの最小値（０として）および最大値（１として）によって正規化された（すなわち、水平軸に沿った正規化）。この分析は、細菌相関ネットワークにおいて他の大部分と高度に相関している細菌（すなわち、ｃＯＴＵ）である「マスター細菌」を見出すことを可能にする。

【図11】強く相関する細菌群（ＳＣＢＧ）の分析。ａ、樹状図の高さの閾値の関数としてのＳＣＢＧの数（図４ｂ）。赤い点線、閾値０．６。ｂ、閾値の関数として最も多くのｃＯＴＵを含むＳＣＢＧ中のｃＯＴＵの数。ｃ、閾値が０．６のとき、ＳＣＢＧにおけるｃＯＴＵ数の分布。ｄ、各ＳＣＢＧにおける試料のＭａ^ｄｉｓｔ１、Ｍａ^ｐｒｏｘ１、Ｍｂ^ｄｉｓｔ、Ｍｂ^ｐｒｏｘ、Ｍｃ^ｄｉｓｔ、およびＭｃ^ｐｒｏｘについてのｃＯＴＵの平均細胞濃度；黒点、ＳＣＢＧ内のすべてのｃＯＴＵは正の相関を示し；紫色および水色の点、全て正の相関を示すｃＯＴＵを示し、異なるサブグループのｃＯＴＵは負の相関を示した。

【図12a】各ＳＣＢＧ内の各ｃＯＴＵと他の間の相対的な相関強度を示すネットワークおよびｒ分布。左、Ｆｏｒｃｅ－ｄｉｒｅｃｔｅｄ ｌａｙｏｕｔ^３９により可視化されたＳＣＢＧのネットワーク。ノード、ｃＯＴＵ；ノード番号、ｃＯＴＵＩＤ；エッジ色、リンクされたｃＯＴＵの間のｒ；ＩＤ色、図１１ｄにおけるドット色と同じ意味。右、ＳＣＢＧにおいて、与えられたｃＯＴＵと他のすべてのｃＯＴＵの間のｒの分布。ｃＯＴＵをまずサブグループ（ＩＤ色）で分け、次いで各サブグループ内の各ｃＯＴＵについて全ての正のｒの平均（青線）で並べられた。

【図12b】同上。

【図12c】同上。

【図12d】同上。

【図12e】同上。

【図12f】同上。

【図13】ＳＣＢＧにおけるｃＯＴＵの門から属への分類。点の色、図１１ｄにおける点色と同じ意味。すべてのＳＣＢＧは複数の属を含み、＞６０％のＳＣＢＧ（１９／３１）は複数の門さえ含んでおり、ＳＣＢＧは分類学とあまり相関していないことを示唆している。一方、≧２ｃＯＴＵを含む検出されたすべての属において、発明者らは複数のＳＣＢＧからｃＯＴＵを見出し、微生物相の細菌ネットワークを理解するためには、属レベルより低い解析、実際にｃＯＴＵレベルが重要であることを示唆した。

【図14】ＳＣＢＧ間の相関。ａ、Ｒ_{ｉｎｎｅｒ}とＲ_{ｉｎｔｅｒ}の分布。ｂ、列、与えられたＳＣＢＧと他のすべてのＳＣＢＧの間のＲ_{ｉｎｔｅｒ}の分布；ＳＣＢＧは、水平軸に沿った分布の平均（青線）によって並べられた。

【図15】シークエンスによって得られた各ｃＯＴＵのカウントに基づく、それぞれＭａ^ｄｉｓｔおよびＭａ^ｐｒｏｘ２サンプルの技術的反復間の比較。Ｍａ、マウス、；ｄｉｓｔ、ｐｒｏｘ、位置；１、２、３、技術的反復、ｒ：Ｐｅａｒｓｏｎ係数。

【図16a】各ｃＯＴＵの細胞濃度の技術的反復間および試料間の比較。Ｍａ^ｄｉｓｔ１－Ｍａ^ｄｉｓｔ３、Ｍａ^ｄｉｓｔ３－Ｍａ^ｐｒｏｘ２、Ｍｂ^ｄｉｓｔ－Ｍｃ^ｐｒｏｘ（赤いアスタリスク記号）の３つの例を図１６ｃに示す。Ｍａ、ＭｂおよびＭｃ、マウス；ｄｉｓｔおよびｐｒｏｘ、位置；１、２、および３、技術的反復。ドット、ｃＯＴＵ；マゼンタ線、ポアソン分布に基づく正規化されたサンプリングノイズの理論的信頼間隔（９９．９％）；水色線、２倍変化；青色ドット、異なる濃度を示すｃＯＴＵ、それぞれ挿入番号、青色および灰色のドットの数。

【図16b】同上。

【図16c】同上。

【図17】ＳＣＢＧのＩＤ。上、図４ｄと同じ樹状図樹であるが、ｃＯＴＵ ＩＤがある。下：赤い正方形の位置、図４ｄに示すヒートマップにおけるＳＣＢＧの位置；青い番号、各ＳＣＢＧに対するＩＤ。

【図18】同じ試料について、プライマーセットＦ１－Ｆｗ／Ｆ１－Ｒｖと３４１Ｆ／８０５Ｒを用いたｄｄＰＣＲ測定の比較。ａ、プライマーＦ１－Ｆｗ／Ｆ１－Ｒｖを用いた盲腸細胞試料のｄｄＰＣＲにより測定した液滴の蛍光強度の分布。ｂ、ａと同じ測定であるが、異なるプライマー３４１Ｆ／８０５Ｒを使用する。ｃ、それぞれｂの蛍光強度分布に当てはめた４種類のガウス分布および４つのガウス分布の混合合計。ｄ、フィッティングしたガウス分布の数の関数としてのフィッティングに基づき計算された陽性の液滴の割合。水色、プライマーＦ１－Ｆｗ／Ｆ１－Ｒｖによって増幅された細胞試料；青、プライマー３４１Ｆ／８０５Ｒによるが、水色と同じ細胞試料；赤、プライマーＦ１－Ｆｗ／Ｆ１－Ｒｖによって増幅された細胞外試料；黒、赤と同じ細胞外試料であるが、プライマー３４１Ｆ／８０５Ｒによって増幅された細胞外試料；エラーバーは、３つの独立したフィッティングの標準偏差（異なるランダムな初期値を有する）。ｅ、同一試料についてプライマーＦ１－Ｆｗ／Ｆ１－Ｒｖを使用する場合とプライマー３４１Ｆ／８０５Ｒを使用する場合のｄｄＰＣＲ測定の比較；４つのガウス分布を用いる適合に基づいて計算された陽性液滴の割合；細胞、細胞試料；ｅｃＤＮＡ、細胞外ＤＮＡ試料；エラーバー、標準偏差、ｎ＝４。

【図19】スパイクインコントロール調製の配列情報を含む概略図。ＳｔｄＴａｒｇｅｔ１、ＳｔｄＴａｒｇｅｔ２、ＲａｎｄｏｍＢａｒ＿ｓｔｄ１、Ｓｔｄ＿Ｒ２、Ｐ２＿ｑＰＣＲ＿Ｒｖ、およびＰ１＿ｑＰＣＲ＿Ｆｗ、合成ＤＮＡオリゴヌクレオチド；ＳｔｄＴａｒｇｅｔ２中の「５Ｐｈｏｓ」、オリゴヌクレオチドの５’末端のリン酸化； Ｉｎｄｅｘは８塩基； 配列中のＮ、Ａ，Ｃ，Ｇ，またはＴ。

【図20】ステップ３．２の論理図。

【図21】ステップ５の論理図。

【図22】Ｓｉｌｖａデータベースに登録されたＶ３－Ｖ４領域における１６Ｓ ｒＲＮＡ遺伝子の長さの分布。プライマー３４１Ｆおよび８０５Ｒにマッチした１６Ｓ ｒＲＮＡ遺伝子のみを用いた（全体の８６．４％）。長さ、３４１Ｒと一致した最初の塩基から８０５Ｒと一致した最後の塩基までの塩基の数。対応する１６Ｓ ｒＲＮＡ遺伝子の全長の総計９９．９４％は、４００～５００の範囲にある。

【図23】ステップ７の論理図。

【図24】ステップ８の論理図。

【図25】ステップ９の論理図。

【図26】Ｓａｎ配列による１塩基の異なるＲｅｐＳｅｑタイプの特性評価。ａ、Ｓａｎ配列に基づくグループ化ＲｅｐＳｅｑタイプの平均カウント（Ｍｏｃｋ－ａ、Ｍｏｃｋ－ｂ、Ｍｏｃｋ－ｃデータに基づく）。塩基の相違、ＲｅｐＳｅｑタイプと各グループで最も近いＳａｎ配列（Ｓａｎ配列をＩＤ群とした）との間の塩基の相違の数。ｂ、各グループで一致したＲｅｐＳｅｑタイプの平均カウントに対するその平均カウントの比に対する１塩基の異なるＲｅｐＳｅｑタイプの最高平均カウント。スカイブルーラベル、グループＩＤ。

【図27】Ｂａｒ配列のｃＯＴＵへのクラスター化。ａ、Ｍｏｃｋ－ｂのデータに基づいたｌｏｇ_１０（Ａ×Ｂ）に対するｌｏｇ_１０（Ｏｖｅｒｌａｐ）。ドット、Ｂａｒ配列の可能なすべてのペア；Ｏｖｅｒｌａｐ、ＡおよびＢは、それぞれＢａｒ配列、ＢＳ＿Ａのみ、およびＢＳ＿Ｂのみを含むＢＣｌｕｓｔｅｒの数である（ＢＳ＿ＡおよびＢＳ＿Ｂはペア内の２つのＢａｒ配列である）。青色の破線、フィッティングの９５％信頼区間。ｂ、Ｍｏｃｋ－ａ、Ｍｏｃｋ－ｂ、Ｍｏｃｋ－ｃのデータに基づいたｌｏｇ_１０（Ａ×Ｂ）＋ＯＤに対するｌｏｇ_１０（Ｏｖｅｒｌａｐ）。ドット、３回のサンプリングで可能なすべてのＢａｒ配列のペア（同一のＢａｒ配列ペアには３つのドットがある）；異なる株、異なる株から同定されたＳａｎ配列と一致したペアのＢａｒ配列；ＪＣＭ／ＡＴＣＣ番号、所定の株から同定されたＳａｎ配列と一致したペアのＢａｒ配列；緑色の線、シミュレーションにより得られたｌｏｇ_１０（Ｏｖｅｒｌａｐ）の分布の９９．９％の片側信頼区間；黄色の線、ｘ＝ｙ；ＯＤ、ａにおけるフィッティングにより推定されたｌｏｇ_１０（Ｄｒｏｐｌｅｔｓ／μ）。このプロットのために異物混入したＢａｒ配列を除いたことに留意されたい（ステップ１４参照）。ｃ、（ｂ）と同様、データはＭ０－ａ、Ｍ０－ｂ、およびＭ０－ｃのデータに基づいていた。各Ｂａｒ配列の名前はＳｉｌｖａデータベースに基づく。異なる名前、１対のＢａｒ配列のマップされた名前が違う；同一名（科）、１対のＢａｒ配列のマップされた名前は同じで、科の名前またはそれより高次の分類の名前しか決められていない；同一名（属）、１対のＢａｒ配列のマップされた名前は同じで、属の名前またはそれより高次の分類名前しか決められていない；不明、１対のＢａｒ配列のうちの片方または両方がデータベースに登録されていない。ｄ、Ｒａｔｉｏ＿Ｐｏｓｉｔｉｖｅの分布（ステップ１２参照）。検体のＭａ^ｄｉｓｔ１－３、Ｍａ^ｐｒｏｘ１－３、Ｍｂ^ｄｉｓｔ、Ｍｂ^ｐｒｏｘ、Ｍｃ^ｄｉｓｔ、Ｍｃ^ｐｒｏｘの結果。

【図28】模擬細胞集団とＭ０の間における検出されたｃＯＴＵの平均カウント（３反復からの）の比較。Ｍ０で検出されなかったｃＯＴＵは示されていない。ＪＣＭ／ＡＴＣＣ、所定の株と一致するｃＯＴＵ。ＣＯＴＵ＜数＞、いずれの設計株とも一致しなかったｃＯＴＵ；Ｉ、ＩＩ、ＩＩＩ、３つのカテゴリー（ステップ１４参照）。

【図29】細菌凝集塊の破壊。ａボルテックス前のＪＣＭ１０１８８の細菌の塊。ｂ、ボルテックス後に１つのドットまたは複数のドットを含むスポットの例。ｃ、ボルテックス後の各菌株及び盲腸試料のスポット当たりの点の数の分布。ｄ、ボルテックス後の各菌株および盲腸試料のスポット当たりの点の数の平均。ｆ、盲腸試料の計２０８個の確認されたスポットのうち複数の点を含むすべてのスポット。黄色の矢印、この一例のみ同一のスポット内に２つの異なる形の点があるように見える。

【図30】顕微鏡画像による細菌数の測定。ａ、同視野の大腸菌（ＤＨ５α）の位相差照明と蛍光照明（ＰＩ）の比較。ｂ、顕微鏡画像による細菌数の概略。ｃ、蛍光照射および位相差で照射し、ＰＩで染色したＡＴＣＣ７００９２６株。バックグラウンド除去の閾値をｅに示した。赤色矢印、位相差照射によっても観察される微小球。ｄ、ｃにおける拡大画像（Ａ－Ｅ）。カラーライン、ＩｍａｇｅＪによる輝度測定に使用された線プロファイル；数、（ｅ）で示される明るいスポット（すなわち、細菌）の番号。ｅ、（ｄ）の線プロファイルに沿って測定された輝度（グレーの値）。破線、背景除去の閾値（ｃ参照）。

【図31】ｅｃＤＮＡと細胞の分離の制御。ａ、孔径０．１μｍ、０．２２μｍ、０．４５μｍのＵｌｔｒａｆｒｅｅ（商標）－ＭＣ遠心フィルターを用いた異なる濾液の比較。フィルター上の残渣、フィルターメンブレンの上に残った試料；フィルター通過液、フィルターメンブレンを通過した液体；存在量、ｄｄＰＣＲで測定した総コピー数。ｂ、ｄｄＰＣＲと顕微鏡画像で測定した、濾過後の細胞とｅｃＤＮＡの存在量を比較したもの。存在量、ｄｄＰＣＲで測定したコピーの総数、または蛍光イメージングで測定した明るいスポットの総数。ｃ、濾過および遠心分離を用いたｅｃＤＮＡと細胞の分離の比較。存在量、ａと同じ。

【図32】盲腸試料の細胞および細胞外ＤＮＡ。ａ、総濃度で正常化した細胞およびｅｃＤＮＡの濃度比率。合計（１００％）は、細胞およびｅｃＤＮＡの濃度の合計として定義した。エラーバー、計算に基づく細胞とｅｃＤＮＡの濃度の標準偏差（ｎ＝５）から計算された伝搬誤差。ｂ、分離した細胞とｅｃＤＮＡの合計濃度の、濾過していない試料の総濃度に対する比。エラーバー、伝搬標準偏差（ｎ＝５）。ｃおよびｄ、各ｃＯＴＵ濃度の細胞およびｅｃＤＮＡの合計濃度と、試料Ｍａ^ｄｉｓｔおよびＭａ^ｐｒｏｘの非濾過試料濃度との比較。赤い点はｅｃＤＮＡが検出されたｃＯＴＵ、黒い点はｅｃＤＮＡが検出されなかったｃＯＴＵ。濾過の繰返し３回はそれぞれ比較した。

【図33】クラスター数（固有のバーコード）のバーコードに設計されたランダムな塩基数への依存性。シークエンスラン１の結果を示す。

【図34】ユニークなバーコードあたりのリード数の平均値について、１０系統のＳａｎ配列とマッチさせたｃＯＴＵの数の依存性。Ｍｏｃｋ－ｂのデータを示した。

【図35】ｃＯＴＵの豊富さ。ａ、各試料について６，０７５細胞をサブサンプリングした。ｂ、３，０００細胞をサブサンプリングした。Ｍａ、ＭｂおよびＭｃ、マウス；ｄｉｓｔおよびｐｒｏｘ、位置；１ａ、２および３、技術的反復；エラーバー、標準誤差。

【図36】ｃＯＴＵの比例的存在量を用いて測定した試料間の微生物相のＢｒａｙ－Ｃｕｒｔｉｓ非類似度。Ｍａ、ＭｂおよびＭｃ、マウス；ｄｉｓｔおよびｐｒｏｘ、位置；１、２、および３、技術的反復。

【図37】同じ分類群からのｃＯＴＵペア（点）のｒの分布（バイオリンプロット）。門から科レベルまでの名前を囲むカラーボックスは、そのオーナーシップを表す。

【図38】ｒに対するｃＯＴＵ間の平均濃度の比の依存性。点、ｃＯＴＵのペア；比、ペアにおいて高い方の濃度を低い方の濃度で割った値；黄色の線、定量的な等高線（１０％間隔）。

【図39】ブレインスライサを用いたマウス大腸試料の細分化のようすを示す。パネルａ～ｆには、マウス大腸試料をブレインスライサ上に載置し（パネルａ）、包埋し（パネルｂ）、凍結させ（パネルｃ）、切断し（パネルｄおよびｅ）、切断（細分化）されたサンプルを得た（パネルｆ）ようすが示されている。パネルｇには、盲腸（Ｃｅｃａｌ側）から肛門（Ａｎａｌ側）に向けての領域分けが図示されている。Ｃエリアに関しては、さらに中心部と周辺部とで細分化された（パネルｇ）。

【図40】各サンプルのバーコード配列の濃度を示す。“－Ｃｅｌｌ”は、細胞非存在下での結果を示し、“＋Ｃｅｌｌ”は、細胞存在下での結果を示す。エラーバーは、標準偏差を示す（ｎ＝４）。

【図41】ｄｄＰＣＲの第三ステージのサイクル数と、液滴の蛍光強度との関係（パネルａ）および全液滴に占める陽性液滴の割合の関係（パネルｂ）を示す。エラーバーは、標準偏差を示す（ｎ＝４）。

【図42】ｄｄＰＣＲの第三ステージの反応時間と、全液滴に占める陽性液滴の割合の関係を示す。エラーバーは、標準偏差を示す（ｎ＝４）。

【発明を実施するための形態】

【0010】

本明細書では、「対象」とは、生物であり、動物および植物である。対象は、例えば、例えば、脊椎動物であり、例えば、哺乳動物、魚類、鳥類、両生類、は虫類、例えば、ヒト、チンパンジー、ゴリラ、オランウータン、サル、マーモセットおよびボノボなどの霊長類、ブタ、ラット、マウス、ウシ、ヒツジ、ヤギ、ウマ、ネコおよびイヌなどの四足動物（例えば、食肉類、偶蹄類、奇蹄類およびげっ歯類）であり得る。

【0011】

本明細書では、「細胞」は、生物の細胞であり、細菌、原生動物、クロミスタ、動物、植物、および菌の細胞であり得る。本明細書では、「単一化された細胞」とは、１つ１つ分離した存在形態を有する細胞を意味する。すなわち、単一化された細胞を含む溶液とは、１以上の細胞を含む溶液であり、それぞれの細胞が１つ１つ分離した存在形態を有する溶液である。単一化された細胞を含む溶液は、好ましくは、含まれる全細胞またはほとんどの細胞が１つ１つ分離された存在形態を示す溶液であるが、単一化された細胞を含む限り、溶液が２つ以上の細胞が接着してなる細胞塊を含んでいてもよい。

【0012】

本明細書では、「細胞集団」は、複数の細胞を含む組成物である。細胞集団は、一般には、複数種類の細胞を含み、各種類は、複数の細胞をそれぞれ含み得る。組成物の形態は、液体または固体であり得る。

【0013】

本明細書では、「微生物叢」とは、微生物の集団である。天然には、様々な微生物叢が存在する。例えば、土、水（海、川、沼、池）、空気、動物の表皮、体毛、口腔、鼻腔、消化管（食道、胃、小腸、大腸、盲腸等）、および生殖器；ならびに、植物の外皮、および根などに微生物叢が存在する。動物における微生物叢は、当該動物の健康状態を反映し、または健康状態に影響し得る。微生物叢は、１０種類以上、２０種類以上、３０種類以上、４０種類以上、５０種類以上、６０種類以上、７０種類以上、８０種類以上、９０種類以上、または１００種類以上の微生物を含み得る。微生物叢は、未知の微生物を含み得る。微生物叢における未知の微生物は、含まれる微生物の種類の１０％以上、２０％以上、３０％以上、または４０％以上であり得る。

【0014】

本明細書では、「細胞バーコード」とは、１つ１つの細胞に対して割り当てられる固有の塩基配列を有する核酸である。１つ１つの細胞が異なる塩基配列（すなわち、細胞に固有の塩基配列）を有する細胞バーコードと紐付けられ得る。従って、細胞バーコードの数は、細胞の数を示し得る。このようにすることで、従来は定量的に測定されてきた細胞数を、定性的に評価可能な塩基配列数に変換して測定することができる。細胞バーコードは、存在する総細胞数に対して十分な種類を用意することができる。

【0015】

本明細書では、「単離」とは、目的物を他から分離することを意味する。単離は、分離した後に、目的物を濃縮または精製することを含み得る。

【0016】

本明細書では、「増幅産物」とは、遺伝子増幅（例えば、ポリメラーゼ連鎖反応（ＰＣＲ））により増幅して得られた核酸をいう。ＰＣＲでは、ＤＮＡの増幅させる部位を挟むように２つのプライマーを設計し、所定の条件下でＤＮＡポリメラーゼと反応させることによって、２つのプライマーで挟まれた部分が増幅される。プライマーは、ＤＮＡの増幅させる部位とハイブリダイズする配列を有する単鎖形態の核酸であり得るが、当該核酸の５’末端には、付加的な塩基配列（例えば、アダプター、サンプルに固有のインデックス配列や制限酵素認識部位等）を連結してもよい。

【0017】

本明細書では、「パラログ」とは、ゲノム上で、遺伝子重複によって生じた２つの遺伝子をいう。本明細書では、「オーソログ」とは、異なる生物に存在する相同な機能を有する遺伝子をいう。

【0018】

本発明によれば、細胞集団の処理方法であって、
（Ａ）単離された細胞集団を含む細胞の分散液から、液滴集団であって、水性の液滴を含み、前記液滴の少なくとも一部のそれぞれが、１つの細胞と１分子の細胞バーコードと
を含む液滴集団を得ること
を含む方法が提供される。

【0019】

本発明によればまた、
液滴集団であって、水性の液滴を含み、前記液滴の少なくとも一部のそれぞれが、１つの細胞と１分子の細胞バーコードとを含む液滴集団が提供される。この態様において、細胞は、単離された細胞集団（例えば、微生物叢）を構成していた細胞であり得る。

【0020】

上記（Ａ）の細胞の分散液は、単離された細胞集団中に含まれる細胞を水溶液中で分散させることによって得られ得る。細胞は、水流、例えば、振とう、ピペッティング等による水によるシアストレスなどによって溶液中に分散させることができる。「分散させる」とは、水溶液中で、複数の細胞を含む細胞塊を乖離させて複数の単一細胞にすること、および、好ましくは単一細胞を水溶液中に浮遊させることを意味する。本発明の方法は、単離された細胞集団中に含まれる細胞を水溶液中で分散させることを含み得る。

【0021】

ある態様では、細胞集団は、微生物叢であり得る。この態様において、微生物叢としては、天然の微生物叢が好ましく用いられ得る。微生物叢としては、例えば、土、水（海、川、沼、池）、空気、動物の表皮、体毛、口腔、鼻腔、消化管（食道、胃、小腸、大腸、盲腸等）、および生殖器；ならびに、植物の外皮、および根に存在する微生物叢を用いることができ、例えば、消化管の微生物叢を用いることができる。微生物叢としては、例えば、口腔内の微生物叢、食道内の微生物叢、胃内の微生物叢、十二指腸内の微生物叢、小腸内（例えば、空腸内または回腸内）の微生物叢、盲腸の微生物叢、大腸内（例えば、上行結腸、横行結腸、下行結腸、Ｓ状結腸、または直腸内）の微生物叢であり得る。天然の微生物叢は、培養せずに分析することが好ましいが、培養に付してから分析することは許容され得る。ある好ましい態様では、微生物叢は、未知の微生物を含む。ある好ましい対象では、未知の微生物の種類は、微生物叢に含まれる微生物の種類の１０％以上、２０％以上、３０％以上、または４０％以上であり得る。ある態様では、細胞集団は、細胞外ＤＮＡを含み得る。細胞外ＤＮＡは、所定の遺伝子を含み得る。細胞外ＤＮＡは、細胞集団を処理する前に、除去してもよい。細胞外ＤＮＡの除去は、後述するようにフィルターろ過または遠心分離により行うことができる。細胞外ＤＮＡは、処理される細胞集団中に含まれていてもよい。

【0022】

細胞集団の単離は、細胞集団を取得することにより行われる。細胞集団の単離は、取得された細胞集団を、細胞以外の１以上の構成要素から分離することをさらに含んでいてもよい。細胞集団を細胞以外の１以上の構成要素から分離することは、フィルターろ過、または遠心分離により行うことができる。フィルターろ過は、例えば、サブμｍの孔径（例えば、０．２２μｍ）を有するフィルターを用いて行うことができ、細胞集団はフィルター上の残渣から回収され得る。

【0023】

本発明では、液滴を作成する前に、単離された細胞集団に含まれる細胞を水溶液中で分散させることができる。ここで、分散させるとは、細胞１つ１つが分離して存在する状態にすることを意味する。分散は、ピペッティングで細胞塊を細胞を破壊しないように崩すことによって達成し得る。水溶液としては、細胞を破壊しない限り特に限定されないが、水、および生理食塩水等を用いることができる。単離された細胞集団は、純水、生理食塩水、および遺伝子増幅用反応溶液等に分散させることができる。

【0024】

ある態様では、液滴は、オイル中で作成することができる。従って、この態様では、（Ａ）で得られる液滴集団は、オイル中に水性の液滴（水滴）を含む。すなわち、（Ａ）で得られる液滴集団は、油中水滴型粒子（油中に分散された水性の液滴集団）であり得る。

【0025】

上記水滴の粒径は、例えば、下限値が１０μｍ～１００μｍの範囲の数値であり得、上限値が５０μ～１０００μｍの数値であり得る。水滴の粒径は、例えば、１０μｍ～１０００μｍ、例えば、２０μｍ～９００μｍ、３０μｍ～８００μｍ、４０μｍ～７００μｍ、５０μｍ～６００μｍ、５０μｍ～５００μｍ、５０μｍ～４００μｍ、５０μｍ～３００μｍ、５０μｍ～２００μｍ、５０μｍ～１５０μｍ、または例えば、約１００μｍであり得る。このような液滴集団は、例えば、マイクロ流体デバイスを用いて当業者であれば適宜作製することができる。このような液滴集団はまた、市販の液滴製造機を用いて作製することができる。市販の液滴製造機としては、例えば、ＢＩＯ－ＲＡＤ社のＱＸ２００ Ｄｒｏｐｌｅｔ Ｇｅｎｅｒａｔｏｒを用いることができる。

【0026】

本発明の細胞集団の処理方法によれば、液滴集団であって、水性の液滴を含み、前記液滴は、１つの細胞と、その細胞に対して固有の１種類の塩基配列を有する細胞バーコード（例えば、ＤＮＡである）の１分子を含む、水性の液滴を含む液滴集団を得ることができる。より具体的には、本発明の細胞集団の処理方法では、例えば、分散された複数の細胞を含む水溶液と、１分子毎に異なる塩基配列を有する細胞バーコードを含む水溶液とを油中で混合することによって１つの細胞、および細胞毎に固有の単一種類の細胞バーコード含む水性の液滴を含む液滴集団が得られる。

【0027】

本発明の細胞集団の処理方法によれば、他の細胞は、その細胞に対して固有の別の１種類の塩基配列を有する細胞バーコードを含む水性の液滴中に含まれる。細胞は、全液滴の５０％以下、４０％以下、３５％以下、３０％以下、２５％以下、または２０％以下（例えば、２０％）に含まれ得る。このようにすることによって、１つの液滴に対して複数の細胞が含まれる確率を低下させることができ、２０％の液滴に細胞が含まれていると仮定したときには、理論上は、細胞を含む液滴の、例えば、９０％以上の液滴が含む細胞の数が１となる。また、細胞バーコードも、全液滴の５０％以下、４０％以下、３５％以下、３０％以下、２５％以下、または２０％以下（例えば、２０％）に含まれ得るようにすることができる。これにより細胞バーコードを含む液滴の、例えば、９０％以上の液滴が含む細胞バーコードの数を１とすることができる。このようにすると、細胞１つと細胞バーコード１分子とを含む液滴が得られ、前記液滴は、全液滴中の１～１０％、２～６％、３～５％、または例えば４％程度となり得る。ある態様では、全液滴に対する細胞を含む液滴の割合を３０％以下（好ましくは約２０％）とし、かつ、全液滴に対する細胞バーコードを含む液滴の割合を３０％以下（好ましくは約２０％）とすることができる。このように、全液滴に対する細胞および細胞バーコードを含む液滴の割合を低下させることで、２つ以上の細胞が１つの液滴に混入する可能性、および２分子以上の細胞バーコードが１つの液滴に混入する可能性を低減するまたは無くすことができる。なお、細胞および細胞バーコードのいずれか一方ならびに両方を含まない液滴の存在は、細胞中の所定の遺伝子と細胞バーコードの連結物を対象としてシークエンスする後の工程には影響しない。

【0028】

得られる液滴集団において、２つ以上の細胞と１つの細胞バーコードを含む液滴の割合は、例えば、０．５％以下、０．４％以下、または０．３％以下であり得、例えば、０．３％～０．５％であり得る。得られる液滴集団において、１つの細胞と２つ以上の細胞バーコードを含む液滴の割合は、例えば、０．５％以下、０．４％以下、または０．３％以下であり得、例えば、０．３％～０．５％であり得る。得られる液滴集団において、２つ以上の細胞と２つ以上の細胞バーコードを含む液滴の割合は、例えば、０．０５％以下、０．０４％以下、または０．０３％以下であり得、例えば、０．０３％～０．０５％であり得る。ここで、２つ以上の細胞または細胞バーコードが含まれる液滴は、少ないほど好ましいが、しかし、そのような液滴が生じることは許容される。

【0029】

水性の液滴は、１つの細胞および１分子の細胞バーコードに加えて、プライマー、および遺伝子増幅用試薬をさらに含んでいてもよい。細胞は、遺伝子増幅反応中に破壊されるので、試薬は、界面活性剤を含んでいる必要はない。また、水性の液滴は、遺伝子増幅反応に適した水溶液（例えば、遺伝子増幅反応溶液）であり得る。

【0030】

オイルは、遺伝子増幅反応（６０℃～１００℃）の環境下で、安定かつ不活性であるオイルであればいずれのオイルでも用いることができる。そのようなオイルとしては、例えば、鉱油（例えば、軽油）、シリコーン油、フッ化油若しくは、他の市販のオイル、またはこれらの組み合わせを挙げることができるが、これらに限定されない。

【0031】

このような条件下で、細胞、細胞バーコード、プライマー、および遺伝子増幅用試薬を含む水溶液から、液滴集団であって、水性の液滴を含み、前記液滴の少なくとも一部のそれぞれが、得られた細胞の１つと１分子の細胞バーコードとを含む液滴集団を得ることができる。より具体的には、細胞、細胞バーコード、プライマー、および遺伝子増幅用試薬を含む遺伝子増幅反応液を作製し、当該溶液から上記の通り液滴集団を得ることができる。

【0032】

本発明によればまた、
細胞集団に含まれる遺伝子配列を決定（または分析）する方法であって、
（Ａ）単離された細胞集団を含む細胞の分散液から、液滴集団であって、水性の液滴を含み、前記液滴の少なくとも一部のそれぞれが、１つの細胞と１分子の細胞バーコードと
を含む液滴集団を得ることと、
（Ｂ）得られたそれぞれの液滴中で、細胞バーコードの増幅産物と所定の遺伝子の増幅産物を得て、さらに、細胞バーコードと所定の遺伝子の全部または一部の塩基配列を含む連結物を得ることと、得られた連結物をシークエンスして所定の遺伝子の塩基配列と細胞バーコードの塩基配列を決定することと
を含む、方法（以下、本発明の配列決定方法という）が提供される。

【0033】

液滴を形成させる際には、細胞バーコードの増幅と細胞内の所定の遺伝子の増幅のためのプライマーセット、ｄＮＴＰ、および耐熱性ＤＮＡポリメラーゼ等のＰＣＲにおいて必要な鋳型以外の構成要素を溶液中（例えば、細胞バーコード溶液中）に混合しておくことで、各液滴にこれらの遺伝子増幅に必要な要素を導入しておくことができる。その後、液滴集団を含む液体をＰＣＲ用のチューブに移し、ＰＣＲによって各液滴内でＤＮＡの増幅反応を誘導することができる。各液滴中での遺伝子増幅により、それぞれの液滴中で、細胞内の所定の遺伝子の増幅産物と、細胞バーコードの増幅産物とを得ることができる。増幅は、例えば、２５サイクル、好ましくは３０サイクル以上の増幅サイクルを含み得る。次いで、それぞれの液滴中で、細胞内の所定の遺伝子の増幅産物と、細胞バーコードの増幅産物とを連結させることができる（例えば、図５参照）。連結は、例えば、細胞バーコードのプライマーの１つとと所定の遺伝子のプライマーの１つを、互いにハイブリダイズ可能な相補的な配列部分を有するように設計することにより、増幅反応（例えば、ＰＣＲ反応）の過程において行うことができる（例えば、図５の配列番号４および５参照）。このようにすることで、１つの細胞に由来する所定の遺伝子の増幅産物の分子それぞれに対して、１種類の細胞バーコードを付与することができる。

【0034】

細胞バーコード（Ｂａｒｃｏｄｅ）は、中央に細胞に固有の塩基配列を有し｛但し、特定番号の塩基配列を配列間で同じ配列としてもよい｝、両端に増幅プライマーがハイブリダイズするための塩基配列を有していることができる。ある態様では、増幅プライマーがハイブリダイズするための塩基配列は、細胞バーコード間での共通配列であり得る。細胞バーコードを増幅するプライマーは、シークエンス用のアダプター配列と上記細胞バーコードの片末端と遺伝子増幅環境下でハイブリダイズすることができる塩基配列を有していることができる。当該細胞バーコードを増幅するプライマーは、サンプルの種別を識別するためのインデックス配列をさらに有していてもよい。細胞バーコードを増幅する別のプライマーは、所定の遺伝子と連結するためのリンカー配列と、上記細胞バーコードの別の片末端と遺伝子増幅環境下でハイブリダイズすることができる塩基配列を有していることができる。

【0035】

所定の遺伝子を増幅するプライマーは、細胞バーコードを増幅するプライマーに含まれるリンカー配列とハイブリダイズするための塩基配列と、所定の遺伝子の増幅する箇所に遺伝子増幅環境下でハイブリダイズするための塩基配列とを有し得る。所定の遺伝子を増幅する別のプライマーは、所定の遺伝子の増幅する箇所に遺伝子増幅環境下でハイブリダイズするための塩基配列とシークエンス用のアダプター配列とを含み得る。当該所定の遺伝子を増幅する別のプライマーは、サンプルの種別を識別するためのサンプルに固有のインデックス配列をさらに有していてもよい。

【0036】

細胞バーコードの増幅産物と所定の遺伝子の増幅産物は、同じリンカー配列を有し、したがって、遺伝子増幅中に、細胞バーコードの増幅産物と所定の遺伝子の増幅産物とを連結したものの増幅産物が得られ得る。

【0037】

シークエンス用のアダプター配列は、シークエンス前のブリッジＰＣＲのための配列を両端に含み得る。シークエンス用のアダプター配列は、シークエンス用プライマーの結合部位を含み得る。シークエンス用のアダプター配列は、サンプルの種別を識別するためのサンプルに固有のインデックス配列を含み得る。ブリッジＰＣＲでは、固相化した２種類のオリゴＤＮＡに対して、シークエンスに供されるＤＮＡであって、両端に上記２種類のオリゴＤＮＡそれぞれにハイブリダイズ可能な配列を有するＤＮＡをハイブリダイズし、その状態で、ＰＣＲによって固相表面上でＤＮＡを増幅させる技術である。

【0038】

したがって、本発明では、水性の液滴を含み、当該液滴は、１つの細胞に由来する所定の遺伝子の増幅産物を含み、前記所定の遺伝子の１分子毎に、当該細胞に固有の１種類の細胞バーコードが連結している、液滴集団もまた、提供される。この液滴集団では、液滴はそれぞれ、異なる１つの細胞に由来する所定の遺伝子と、当該細胞に固有の１種類の細胞バーコード｛すなわち、液滴毎に異なる細胞バーコードが含まれる｝を含む。

【0039】

１つの細胞に由来する所定の遺伝子の１分子毎に、当該細胞に固有の１種類の細胞バーコードが連結した連結物は、上記のように、シークエンス用のアダプター配列、細胞バーコード配列、リンカー配列、所定の遺伝子の全部または一部の塩基配列、およびシークエンス用のアダプター配列を、この順番で含み得る。この連結物は、サンプルに固有の塩基配列を有するインデックス配列をさらに含んでいてもよい。インデックス配列は、シークエンス用のアダプター配列、細胞バーコード配列、リンカー配列、所定の遺伝子の全部または一部の塩基配列、およびシークエンス用のアダプター配列のいずれか２つの間に含まれ得る。インデックス配列は、代わりに、または追加で、シークエンス用のアダプター配列中に含まれていてもよい。

【0040】

本発明では、１つの細胞に由来する所定の遺伝子の増幅産物の分子それぞれと細胞毎に固有の１種類の細胞バーコードの増幅産物との連結物を作製することができる。ここで、所定の遺伝子は、１種類であることが好ましいが、１種類とは限らず、複数種類である場合もあり得る。細胞バーコードは細胞毎に１種類であることが好ましい。

【0041】

本発明では、ある連結物に関して、決定された所定の遺伝子の塩基配列と細胞バーコードの塩基配列は、紐付けて管理されている。この紐付けに基づいて、同一の細胞バーコードが連結された所定の遺伝子は、同一細胞に由来するものと推定され得る。したがって、本発明の配列決定方法は、それぞれの連結物に関して、決定された所定の遺伝子の塩基配列と細胞バーコードの塩基配列を含む塩基配列の組合せを得ることをさらに含み得る。

【0042】

また、本発明の配列決定方法は、同一の細胞バーコードが連結された所定の遺伝子が同一の細胞に由来すると推定することをさらに含み得る。

【0043】

本発明のある態様では、所定の遺伝子は、微生物の内在性遺伝子であり、好ましくは、進化上、広く様々な種が共有する遺伝子、例えば、ハウスキーピング遺伝子であり得る。ハウスキーピング遺伝子は、エネルギー代謝や細胞機能に必須な遺伝子であり、あらゆる細胞で発現するまたは発現するであろう遺伝子である。ハウスキーピング遺伝子としては、特に限定されないが、例えば、リボソーマルＲＮＡ（ｒＲＮＡ、例えば、１６Ｓ ｒＲＮＡ、および２３Ｓ ｒＲＮＡ）、１６Ｓ ｒＲＮＡと２３Ｓ ｒＲＮＡとの間に存在するｒｉｂｏｓｏｍａｌ ｉｎｔｅｒｇｅｎｉｃ ｔｒａｎｓｃｒｉｂｅｄ ｓｐａｃｅｒｓ（ＩＴＳ）、推定ＡＢＣトランスポート（ａｂｃＺ）、アデニレートキナーゼ（ａｄｋ）、シキメートデヒドロゲナーゼ（ａｒｏＥ）、グルコース－６－ホスフェートデヒドロゲナーゼ（ｇｄｈ）、単一機能ペプチドグリカントランスグリコシラーゼ（ｍｔｇ）、推定デヒドロゲナーゼサブユニット（ｐｄｈＣ）、ホスホグルコムターゼ（ｐｇｍ）、ｒｅｇｕｌａｔｏｒ ｏｆ ｐｉｌｉｎ ｓｙｎｔｈｅｓｉｓ（ｐｉｌＡ）、プロリンイミノペプチダーゼ（ｐｉｐ）、ポリホスフェートキナーゼ（ｐｐｋ）、および３－ホスホセリンアミノトランスフェラーゼ（ｓｅｒＣ）が挙げられる（Ｍａｉｄｅｎ ｅｔ ａｌ．， ＰＮＡＳ， Ｖｏｌ．９５， ３１４０－３１４５， １９９８参照）。これらの遺伝子の配列は、微生物叢の解析において用いられ得る。また、真菌の解析においては、１８Ｓ ｒＲＮＡも用いられ得る。所定の遺伝子が、２種類以上の遺伝子である場合、それぞれの遺伝子が細胞バーコードと連結されるように適切なプライマーおよび反応条件下で増幅反応が行われる。本発明の方法では、所定の遺伝子の塩基配列に基づいて細胞集団を分析するため、できるだけ多くの細胞が有する遺伝子を所定の遺伝子として用いることが有利である。本発明のある態様では、所定の遺伝子は、１６Ｓ ｒＲＮＡをコードする遺伝子であり得る。所定の遺伝子の塩基配列は、当該遺伝子の全長、またはその部分配列であり得る。例えば、１６Ｓ ｒＲＮＡの場合には、シークエンスするのは、全長でなくてもよく、その一部としてもよい。１６Ｓ ｒＲＮＡの一部は、Ｖ３領域およびＶ４領域とすることができる。

【0044】

本発明において、所定の遺伝子は１種類の遺伝子（または相同な遺伝子群）のみを用いれば足り、２種類以上の異なる遺伝子（または互いに非相同な２つ以上の遺伝子群）を用いる必要は無い。但し、所定の遺伝子は２種類以上の異なる遺伝子（または互いに非相同な２つ以上の遺伝子群）であってもよい。

【0045】

本発明の配列決定方法では、シークエンスは、液滴を破壊し、すべての液滴に含まれる溶液を混合して行うことができる。本発明の配列決定方法では、シークエンスは、当業者に周知の方法を用いて実施できる。例えば、シークエンスは、次世代シークエンサー（例えば、イルミナ社のＭｉＳｅｑやＨｉＳｅｑ）を用いて並列的に行うことができる。このように並列に解読するシークエンサーを用いることで、数万～数億の遺伝子断片を迅速に分析することができる。この場合には、当業者は、シークエンスに必要であれば、シークエンス用のアダプターを連結物に付加することができ、この工程は当業者に周知である。

【0046】

本発明の配列決定方法は、シークエンス前に、溶液のＤＮＡを回収することをさらに含んでいてもよい。ＤＮＡの回収は、液滴それぞれに分離して含まれる水相を回収することにより行われ得る。例えば、ＤＮＡの回収は、得られた液滴集団を有機溶媒（例えば、クロロホルムなど）、および、好ましくはさらに水溶液（例えば、緩衝液、例えば、二価金属イオンキレーター（例えば、Ｃａ^２＋のキレーターおよびＭｇ^２＋のキレーター、例えば、エチレンジアミン四酢酸（ＥＤＴＡ））を含むＴｒｉｓ緩衝液、すなわちＴｒｉｓ－ＥＤＴＡ緩衝液またはＴＥ溶液）を添加して、十分に攪拌して、水相と有機相とを分離させ、水相を回収することによって行われ得る。これにより、油中水滴型粒子において、液滴のコンパートメントそれぞれに離散的に存在していた目的のＤＮＡ（すなわち、連結物）を水溶液中に回収することができる。このようにして得られた水溶液では、含まれる液滴に由来するすべての連結物が溶液中（オイルによる区画を有しない溶液である）で混在している（すなわち、液滴のコンパートメントそれぞれに離散的に存在していた連結物が１つの溶液コンパートメント中に存在する状態になる）。上記のようにシークエンスでは、多数の遺伝子断片の塩基配列を並列的に解読できるので、多数のＤＮＡが混在した溶液は、シークエンスに適する。

【0047】

本発明の配列決定方法はまた、シークエンス前に、ＤＮＡを精製することをさらに含んでいてもよい。ＤＮＡの精製は、上記回収工程により得られた水溶液をゲル濾過することにより行われ得る。ゲル濾過としては、ＤＮＡ増幅産物とそれ以外の溶液中の構成要素（例えば、未連結のバーコード増幅物、増幅に用いられなかったプライマー、その他）を分離するために通常用いられる手法によりゲル濾過カラム等を用いて行うことができる。ゲル濾過カラムとしては、例えば、ＤＮＡ精製用のゲル濾過カラムを用いることができる。また、本発明の配列決定方法は、カルボキシル基を表出するカラムまたはビーズで溶液中に含まれるＤＮＡを精製することをさらに含んでいてもよい。カルボキシル基を表出するカラムまたはビーズには脱水和したＤＮＡが塩を介して特異的に吸着することができ、その後、水和させることによりＤＮＡをカラムから離脱させることができる。カルボキシル基を表出するビーズとしては、例えば、Ａｇｅｎｃｏｕｒｔ ＡＭＰｕｒｅ ＸＰ（ベックマンコールター）などを使用することができる。

【0048】

さらには、また、ＤＮＡ増幅の工程において、ＤＮＡの増幅反応を、タグを有するタグ化プライマー（例えば、ビオチン化プライマー）を用いて行った場合には、ＤＮＡ増幅産物はタグ（例えば、ビオチン）が結合している。そのようなタグ付きＤＮＡ増幅産物は、タグに結合する分子（例えば、アビジン、ストレプトアビジン、およびニュートラビジン等のタグ結合分子）を連結したカラムやビーズにより濃縮または除去することができる。特に、本発明の配列決定方法は、細胞バーコードと所定の遺伝子とが連結した産物を得る場合に、連結し損ねたものを除去することを好ましく含み得る。すなわち、細胞バーコードを増幅するプライマーの１つと、所定の遺伝子を増幅するプライマーの１つを、タグを有するプライマーとし、かつ、相補的な配列を有するように設計することができる。すなわち、相補的な配列を有するように設計されたプライマーそれぞれについてのみタグを付与することができる。このようにすると、図５に示されるように、細胞バーコードの増幅産物と所定の遺伝子の増幅産物は、タグ化したプライマー部分に対応する領域において互いに連結し得る。得られた連結物をさらに遺伝子増幅により増幅すると、連結物の増幅物にはタグは含まれないのに対して、連結し損ねた増幅産物はその末端にプライマー由来のタグを有する。このように、ある態様では、前記（Ｂ）において、細胞バーコードを増幅する２つのプライマーのうちの１つであって、所定の遺伝子を増幅する２つのプライマーの１つ（当該１つはタグ分子を有する）と相補的な配列を有するプライマーは、タグ分子を有する。タグ分子は、細胞バーコードと所定の遺伝子の増幅中に、細胞バーコードと所定の遺伝子との連結物からは失われ、タグ分子は、連結し損ねた増幅物においてのみ残ることとなる。したがって、連結し損ねた増幅産物は、タグ結合分子を結合したカラムやビーズを用いてアフィニティーにより除去することができ、これにより、連結物の増幅物をより純度高く精製することができる。
したがって、本発明の配列決定方法は、
前記（Ｂ）において、細胞バーコードの増幅産物は、第一のプライマーに由来する第一の領域を有し、所定の遺伝子の増幅産物は、第二プライマーに由来する第二の領域を有し、第一の領域と第二の領域は、互いにハイブリダイズ可能な相補的な配列部分を有し、前記第一のプライマーおよび第二のプライマーはそれぞれ、１以上のタグ分子を連結しており、当該タグ分子は、前記連結物には含まれず、かつ、
前記（Ｂ）において、水溶液中に回収された連結物から、タグ分子を有する増幅産物を当該タグ分子に親和性を有する分子を担持したカラムまたはビーズを用いて除去することをさらに含んでいてもよい。これにより、所望の連結物から、タグ分子を有する連結し損ねた増幅産物を分離することができる。

【0049】

本発明の配列決定方法では、シークエンス品質の低い塩基配列の領域を削除することを含んでいてもよい。シークエンスの品質は、例えば、Ｐｈｒｅｄアルゴリズムに基づく品質スコア（例えば、ｐｈｒｅｄ品質スコア、例えば、Ｑスコア（Ｑ＝－１０ｌｏｇ_１０（ｅ）｛ここで、ｅは、ベースコールが誤っている確率の推定値である｝）によって行うことができる（Ｅｗｉｎｇ ｅｔ ａｌ．， Ｇｅｎｏｍｅ Ｒｅｓ．， ８（３）： １７５－１８５， １９９８、およびＥｗｉｎｇ ａｎｄ Ｇｒｅｅｎ， Ｇｅｎｏｍｅ Ｒｅｓ．， ８（３）： １８６－１９４， １９９８参照）。シークエンスの解読エラーを低減させるために当業者により広く行われているように、品質スコアが一定の閾値以下である配列は、解析から除外することができる。例えば、Ｑスコアが、２０以下、１５以下、または１０以下である配列を解析から除外することができる。

【0050】

本発明の配列決定方法は、
（Ｃ－１）細胞バーコードの塩基配列に基づいて、決定された塩基配列をクラスタリングして、複数の第一のクラスターを得ること
をさらに含んでいてもよい。

【0051】

上記（Ｃ－１）において「決定された塩基配列」とは、決定された所定の遺伝子の塩基配列と細胞バーコードの塩基配列を含む塩基配列の組合せであり得る。

【0052】

決定された細胞バーコードの塩基配列に基づいて、決定された塩基配列をクラスタリングすることは、細胞バーコードの塩基配列が完全に同一の配列であるか否かでクラスタリングすることのみならず、多少の相違を有する配列を同じクラスターにクラスタリングすることを含み得る。多少の相違を有する配列を同じクラスターにクラスタリングする理由は、実験的には、細胞バーコードの増幅反応やシークエンスにおいてエラーが生じ、解読された塩基配列が本来の塩基配列と異なる配列になり得るためである。しかし、増幅反応やシークエンスにおいて生じるエラーは、経験的によく知られており、それに対応するエラーによって同一配列を異なる配列として区別しないようにするためには、多少の相違を有する配列を同じクラスターにクラスタリングすることが有効である。

【0053】

例えば、決定された細胞バーコードの塩基配列が、完全に同一であるという基準（距離０）でクラスタリングをすると、増幅エラーやシークエンスエラーが存在しない場合には、１つの細胞に由来する塩基配列は、１つのクラスターに正しくクラスタリングされることとなる。したがって、このようなケースについては、問題がない。これに対して、決定された細胞バーコードの塩基配列が、完全に同一であるという基準（距離０）でクラスタリングをすると、増幅エラーやシークエンスエラーが存在する場合には、１つの細胞に由来する塩基配列が、異なる細胞に由来するものとして、２つ以上のクラスターに誤ってクラスタリングされ得る。

【0054】

理論的には、ｎ個の塩基の付加、削除、欠失、または挿入（特にインデル）を有する配列も同一のクラスターにクラスタリングするという基準（距離ｎであり、ｎは１～５の自然数であり得る）では、増幅エラーやシークエンスエラーによってｎ個までの塩基の付加、削除、欠失、または挿入（特にインデル）が発生した場合にも、１つの細胞に由来する塩基配列は、１つのクラスターに正しくクラスタリングされることとなる。ここで、当業者であれば、増幅反応におけるエラー率やシークエンスエラー率に基づいて、ｎを適宜設定することができる。ｎを大きく設定する場合には、細胞バーコードが細胞毎にｎ塩基よりも必ず多く異なるようにデザインすることができる。本発明のある態様では、ｎは、１とすることができる。本発明の別の態様では、ｎは、２とすることができる、本発明のさらなる別の態様では、ｎは、３とすることができる。決定された細胞バーコードの塩基配列が、細胞毎に大きく異なる配列となるように設計することによって、ｎ個までの塩基の付加、削除、欠失、または挿入（特にインデル）が発生した場合にも、当該エラーを有する塩基配列が、いずれの細胞バーコードに由来するかを決定することができる。クラスタリングにはこのような実験的エラーによる影響を低減する効果が期待できる。クラスタリングに関しては、引用することによりその全体が本明細書に組込まれるＷＯ２０１８／２３５９３８を参照して実施することができる。

【0055】

細胞バーコードは、各細胞に対して固有の配列であるために、同一の細胞バーコードを含む連結物は、理論的には、同一細胞に由来する所定の遺伝子としか連結していないはずである。したがって、同一の細胞バーコードの塩基配列に基づいて、決定された塩基配列（細胞バーコードと所定の遺伝子の増幅産物を含む）をクラスタリングすることで、同一細胞に由来する所定の遺伝子が決定できる。所定の遺伝子が、細胞内に１つのみ存在する場合は、理論上、上記（Ｃ－１）で得られる第一のクラスターには、所定の遺伝子に関して１つの配列しか検出されない。他方で、所定の遺伝子が、細胞内に複数存在する場合には、理論上、上記（Ｃ－１）で得られる第一のクラスターは、所定の遺伝子に関して２つ以上の配列（パラログ）を含み得る。そのため、上記（Ｃ－１）をさらに含む本発明の細胞集団の分析方法では、細胞集団中に、所定の遺伝子の重複（コピーまたはパラログ等）を有する細胞の存在を検出することができる。

【0056】

また、上記において、本発明の方法では、細胞数は、理論上は、細胞バーコードの種類の数または細胞バーコードの塩基配列に基づいて得られたクラスターの数に等しいと推定できる。したがって、１細胞中での所定の遺伝子の重複は、算出される細胞数の正確性に影響を及ぼさないという利点がある。

【0057】

したがって、本発明の細胞集団に含まれる遺伝子配列の決定方法は、
（Ｄ－１）得られた第一のクラスターの数から細胞集団に含まれる細胞の数または特定の所定の遺伝子を有する細胞の数を推定すること
をさらに含んでいてもよい。

【0058】

また、上記においては、細胞バーコードの塩基配列に基づいて、決定された塩基配列をクラスタリングして、複数の第一のクラスターを得た。これに対して、以下の実施態様では、本発明の細胞集団を分析する方法は、決定された所定の遺伝子の塩基配列に基づいて決定された塩基配列をクラスタリングすることを含むことができる。

【0059】

すなわち、本発明の細胞集団を分析する方法は、
（Ｃ－２）決定された所定の遺伝子の塩基配列に基づいて、決定された塩基配列をクラスタリングして、複数の第二のクラスターを得ること
をさらに含んでいてもよい。

【0060】

所定の遺伝子の塩基配列に基づいて、決定された塩基配列をクラスタリングして、複数の第二のクラスターを得ることは、完全に同一の配列であるか否かでクラスタリングすることのみならず、多少の相違を有する配列を同じクラスターにクラスタリングすることを含み得る。多少の相違を有する配列を同じクラスターにクラスタリングする理由は、実験的には、細胞バーコードの増幅反応やシークエンスにおいて配列にエラーが生じ得るためである。

【0061】

例えば、決定された所定の遺伝子の塩基配列が、完全に同一であるという基準（距離０）でクラスタリングをすると、増幅エラーやシークエンスエラーが存在しない場合には、クラスタリングにより得られた各クラスターは、正確に１種類の遺伝子に由来することとなる。一方で、理論的には、ｎ個の塩基の付加、削除、欠失、または挿入（特にインデル）を有する配列も同一のクラスターにクラスタリングするという基準（距離ｎであり、ｎは１～５の自然数であり得る）では、増幅エラーやシークエンスエラーによってｎ個までの塩基の付加、削除、欠失、または挿入（特にインデル）が発生した場合にも、１種類の遺伝子に由来するものが同一クラスターにクラスタリングされることになる。そして、本発明のある態様では、ｎは、１とすることができる。本発明の別の態様では、ｎは、２とすることができる、本発明のさらなる別の態様では、ｎは、３とすることができる。ここで、当業者であれば、ｎを適宜設定することができる。得られたクラスターの数は、所定の遺伝子の種類の数に対応する。

【0062】

所定の遺伝子の配列が、すべての微生物において把握されているわけではない。しかし、本発明の配列決定方法では、細胞集団が、未知の微生物を含んでいてもよい。当該未知の微生物が、他の微生物と区別できる塩基配列を有する所定の遺伝子を有する限り、当該未知の微生物を既知の微生物と異なるものとして処理できるためである。

【0063】

ところで、未知の微生物において、所定の遺伝子の塩基配列が、既知の所定の遺伝子の配列と距離ｎ以下の相違しか有しない場合には、上記方法では、未知の遺伝子と既知の遺伝子が本来的に異なる塩基配列を有する場合であっても、同じクラスターにクラスタリングされることとなり、同一の遺伝子に由来すると推定されることとなる可能性がある。

【0064】

したがって、上記（Ｃ－２）は、さらなる工程：
（Ｃ－２α）あるクラスターにおいて、所定の遺伝子に関して相違する塩基配列が含まれている場合に、当該相違する塩基配列の１つの位置において、最も豊富な塩基を決定することと、二番目に豊富な塩基を決定することと、当該位置において、最も豊富な塩基を有する塩基配列の数（すなわち、リード数）に対する二番目に豊富な塩基を有する塩基配列の数（すなわち、リード数）の比（Ｒａｔｉｏ２ｎｄ）が所定の値以上である場合には、最も豊富な塩基を有する塩基配列と、二番目に豊富な塩基を有する塩基配列とを別のクラスターにクラスタリングすること
をさらに含み得る。これによって、同一クラスターに分類された塩基配列のうち、本来的に異なる遺伝子に由来するものを異なるものとして処理することができ、これにより、（ｃ－２）の工程によって、異なる遺伝子が同一と評価される頻度を低減することができる。

【0065】

工程（Ｃ－２α）は、すべての塩基配列の相違について、Ｒａｔｉｏ２ｎｄが所定の値未満となるまで続けられ得る。所定の値は、例えば、０．６以上、０．６５以上、０．７以上、０．７５以上、または０．８以上の数字であり得る。これは、本当に存在する塩基配列であれば、複数の細胞に含まれているはずであり、一定の割合を占めるように検出されるであろうからである。一方で、エラーはその頻度が低いので、この評価によってエラーと本来的に存在する配列とを区別することができる。

【0066】

工程（Ｃ－２α）では、所定の遺伝子の塩基配列についての品質スコアにより、上記リード数に重み付けを行ってもよい。品質スコアは、例えば、Ｐｈｒｅｄアルゴリズムに基づいて決定され得るスコア、例えば、ｐｈｒｅｄ品質スコア、または例えば、Ｑスコアであり得る。品質スコアが所定の値未満である場合には、低く（例えば、０と）重み付けし、品質スコアが所定の値以上である場合には、高く（例えば、スコアの数値に応じて）重み付けしてもよい。実施例のステップ３．２に記載される通りである。

【0067】

この工程によって、最も豊富な塩基を有する塩基配列を当該クラスターにおける「代表的な塩基配列」（ＲｅｐＳｅｑ）とする。

【0068】

異なるＲｅｐＳｅｑを対比し、塩基のずれ（シフト）が見出された場合（すなわち、塩基配列をずらすと２つの塩基配列が一致する場合）には、より多くの第一のクラスターにおいて見出されたＲｅｐＳｅｑを母親とし、より少ない方をシフトとし、シフトした塩基配列を削除して、母親の塩基配列を有すると推定することができる。この際に、シフトした塩基配列のカウント（リード数）を母親であるＲｅｐＳｅｑのリード数に加えることができる。実施例のステップ５に記載される通りである。

【0069】

工程（Ｃ－２α）では、単一のリードでしか検出されなかった塩基配列は、エラーとして除外することをさらに含んでいてもよい。

【0070】

複数の第一のクラスターから同じ配列が検出される場合、その配列は、真に存在していた配列である可能性がある。そのため、工程（Ｃ－１）と工程（Ｃ－２）とは組み合わせて実施することにより、さらに塩基配列の決定精度が高まる。また、工程（Ｃ－１）と工程（Ｃ－２）とを組み合わせると、１つの細胞に複数の所定の遺伝子が存在したかどうかを決定することができる。

【0071】

これによって、本発明では、
（Ｄ－２）得られた第二のクラスターの数から細胞集団に含まれる細胞の種類の数（何種類の細胞が細胞集団に含まれるか）を推定すること
ができる。

【0072】

したがって、本発明の細胞集団に含まれる遺伝子配列を決定する方法は、
（Ｃ－３）決定された細胞バーコードの塩基配列に基づいて、決定された塩基配列をクラスタリングして、複数の第一のクラスターを得ることと、決定された所定の遺伝子の塩基配列に基づいて、決定された塩基配列をクラスタリングして、複数の第二のクラスターを得ることと
をさらに含んでいてもよい。

【0073】

ここでのクラスタリングの工程の詳細およびその効果は、上記（Ｃ－１）および（Ｃ－２）で説明した通りである。上記（Ｃ－３）では、第一のクラスターそれぞれに対して、第二のクラスターを形成させてもよく、第二のクラスターそれぞれに対して第一のクラスターを形成させてもよい。

【0074】

本発明の細胞集団に含まれる遺伝子配列を決定する方法は、
（Ｄ－３）得られた細胞バーコードの塩基配列と所定の遺伝子の塩基配列の組合せの情報に基づいて、少なくとも１つのある第二のクラスターに分類された所定の遺伝子の塩基配列と連結している細胞バーコードの塩基配列から当該所定の遺伝子の塩基配列が分類された第一のクラスターを決定し、当該細胞バーコードが分類された第一のクラスターの数から、当該第二のクラスターに分類された細胞の数を推定すること
をさらに含んでいてもよい。

【0075】

ここで、所定の遺伝子の塩基配列に関して、距離ｎの２つの塩基配列（例えば、配列の中心部において１つの喪失または欠失の相違（すなわち、１－インデル）を有する２つの塩基配列）について、より多くの第一のクラスターに分類された塩基配列を母親とし（すなわち、より多くの数の細胞において検出される塩基配列を母親とし）、より少ない方を１－インデルとする。母親のリード数の方が１－インデルより多い第一のクラスターの数（Ｎｏ_{ｍｏｔｈｅｒ}）と、母親のリード数が１－インデルのリード数よりも少ない第一のクラスターの数（Ｎｏ_{１－インデル}）とを比較して、Ｎｏ_{ｍｏｔｈｅｒ}がＮｏ_{１－インデル}よりも大きい場合には、その母親と１－インデルのペアを残すことができる。さらに、母親と１－インデルの両方を含む第一のクラスターの数に対する、Ｎｏ_{１－インデル}の比が、所定の値（例えば、（Ｎｏ_{１－インデル}－３）／Ｎｏ_{１－インデル}）よりも小さい場合に、母親と１－インデルのペアを残すことができる。残った母親と１－インデルのペアに対して、１－インデルを削除し、１－インデルのリード数を母親のリード数に足すことができる。また、同じ１－インデルに対して母親が２通り存在する場合には、より多くの第一のクラスターにおいて見出された母親に対してそのリード数を足すことができる。また、もし、母親なしで１－インデルのみが検出された第一のクラスターが存在した場合には、そのクラスターにおいては、１－インデルのリード数を母親のリード数とすることができる。実施例のステップ７に記載される通りである。

【0076】

さらに、遺伝子増幅の過程において、ある増幅産物が他の増幅産物と連結したキメラ分子が生じることが問題となり得る。実施例において、本発明の方法では、キメラ分子の生成割合が極めて低いことが明らかであるが、本発明の方法は、当該キメラ分子を特定することをさらに含み得る。キメラ分子は以下のように特定することができる。例えば、キメラ分子を含む第一のクラスター数とキメラ分子（Ｔｏｔａｌ＿Ｎ）のみを含み、親を含まない第一のクラスター数（Ｎ＿ｄ）との比（Ｎ＿ｄ／Ｔｏｔａｌ＿Ｎ）の比が１未満の一定値以下である場合には、このキメラ分子は、エラーによって生じたとしてＲｅｐＳｅｑから除外することができる。実施例のステップ８に記載される通りである。

【0077】

本発明の方法は、上記（Ｃ－２）に加えて、
細胞ベースの操作上分類単位（ｃＯＴＵ）を作成することをさらに含んでいてもよい。細胞集団に含まれる微生物数及び種類が不明であることが多く、さらには、未知の微生物が存在する場合には、データベースに登録された既知の遺伝子配列情報からのみでは、細胞集団の遺伝子配列の分析は不十分となる。特に、所定の遺伝子の塩基配列ベースで操作上部類単位（ＯＴＵ）を形成させると、所定の遺伝子についてある微生物種でｎ個の重複を有する場合には、当該微生物種の数が本来のｎ倍あるものとしてカウントされて誤差を生じることになる。また、２つの異なる微生物種において、一方は塩基配列ＡとＢとを有し、他方はＡとＣとを有する場合に、塩基配列ベースで操作上部類単位（ＯＴＵ）を形成させると、ＯＴＵは、ＡとＢとＣそれぞれに対応して３つ形成され、Ａのカウントの分だけ細胞数に誤差を生じることになる。そこで、（Ｃ－４）では、遺伝子重複を有する細胞が細胞集団に含まれるときのカウントの上記誤差を低減するために、ＲｅｐＳｅｑの情報からｃＯＴＵを作成する。なお、ｃＯＴＵは、理論的には、所定の遺伝子の塩基配列で分類できる微生物の分類単位であり、これまで上位の分類群でしか分類できなかった微生物をさらに詳細に分類する技術的手段である。これは、特に詳細な分類がなされていない微生物や未同定の微生物を含む細胞集団の分析において有用である。分類ができれば、これを元にして細胞集団間での相違を比較することができ、有利である。

【0078】

ｃＯＴＵの作成は以下のように行うことができる。すなわち、従来法同様に、一つの第二のクラスターを一つのｃＯＴＵとみなすことができる。しかし、本発明では、同じ細胞に２つ以上の第二のクラスターが含まれることを考慮して、同一の細胞バーコードと連結された複数の第二のクラスターは１つのｃＯＴＵに分類することをさらに含むことができる。

【0079】

すなわち、本発明の方法は、例えば、上記（Ｃ－３）に加えて、
（Ｃ－４）同一の第一のクラスターに分類された配列が異なる第二のクラスターに分類される場合、当該第二のクラスターを同一の細胞ベースの操作上分類単位（ｃＯＴＵ、すなわち同一の細胞分類）に分類すること
をさらに含んでいてもよい。

【0080】

このｃＯＴＵの作成においては、実験的エラー（例えば、１つの液滴中に２つの細胞が含まれて解析され、これにより一つの第一のクラスターに２つの細胞由来の所定の遺伝子の塩基配列が検出されること）を除外することをさらに含んでいてもよい。

【0081】

ある細胞バーコード配列と連結した所定の遺伝子の塩基配列が２つ存在する場合に、１つの液滴に２つの細胞が混入する確率は、理論上、ポアソン分布にしたがう。上記エラータイプＡは、液滴作成時の細胞の濃度に依存するエラーであると考えられるために、液滴作成時の濃度を薄めること（実施例では２０％の液滴に細胞が含まれるような濃度を用いた）、によって、このエラーの頻度を低下させることができると考えられる。また、２つの塩基配列が異なる細胞に存在していたときに操作中に１つの液滴内に含まれる確率は、理論上、ポアソン分布に従う。

【0082】

１つの細胞バーコードが付与されたＲｅｐＳｅｑ（ＲｅｐＳｅｑは、上記のステップで様々なエラー除去後の配列としてもよく、それが好ましい）が複数存在する場合には、それらをすべてピックアップする。２つのＲｅｐＳｅｑが含まれる液滴の数は（Ｏｖｅｒｌａｐ）と表す。異なる細胞に由来する２つのＲｅｐＳｅｑが１つの液滴に含まれる確率（Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）は、（Ａ×Ｂ×μ）／液滴総数と表される｛ここで、細胞総数は、細胞バーコードを含む液滴の総数であり、Ａは一方のＲｅｐＳｅｑが含まれる液滴の数であり、Ｂは他方のＲｅｐＳｅｑが含まれる液滴の数であり、μは、ＰＣＲ増幅効率、シークエンス深度効果などを含み得る液滴における検出効率のための統合パラメータである｝。ここで、
式：（Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）＝（Ａ×Ｂ×μ）／液滴総数は、
ｌｏｇ_１０（Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）＝ｌｏｇ_１０｛（Ａ×Ｂ×μ）／液滴総数｝
に変換することができる。さらに上記式は、
ｌｏｇ_１０（Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）
＝ｌｏｇ_１０（Ａ×Ｂ）－ｌｏｇ_１０（液滴総数／μ）
変換することができる。ここで、ＡおよびＢは実験的に測定でき、ｌｏｇ_１０（液滴総数／μ）は実験毎に一定の定数とすることができる。従って、ｌｏｇ_１０（液滴総数／μ）を定数ＯＤとすると、上記式は、
ｌｏｇ_１０（Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）＝ｌｏｇ_１０（Ａ×Ｂ）－ＯＤ
に変換できる。これは、ｙ＝ｘ－ＯＤで直線近似することができる。ＡとＢについて様々な整数を想定してｌｏｇ_１０（Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）を算出することができる。現実のｌｏｇ_１０（Ｏｖｅｒｌａｐ）の値が、計算されるｌｏｇ_１０（Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）の信頼区間の外である場合、２つの塩基配列は１つの細胞内に含まれていたと推測することができる。また、ｌｏｇ_１０（Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）の信頼区間の内部である場合、２つの塩基配列はそれぞれ異なる細胞に含まれていたと推測することができる。信頼区間としては、例えば、片側信頼区間（例えば、９５％以上、９８％以上、９９％以上、または９９．９％もしくはそれ以上の信頼区間とすることができる）を用いることができる。このようにして、統計学的にポワソン分布では説明ができない場合に、２つの塩基配列は１つの細胞内に含まれていたと推定することができる。あるいは、統計学的にポワソン分布で説明できるときには、２つの塩基配列は異なる細胞に存在したと推定することができる。

【0083】

また、理論的には、同じ微生物中のＲｅｐＳｅｑの結果は、異なるサンプルにおいても同じであると考えられる。したがって、異なるサンプルにおいても、再現された場合には、複数の異なる細胞集団サンプルを測定し、２つのＲｅｐＳｅｑを含むサンプル数に対する、ｌｏｇ_１０（Ｏｖｅｒｌａｐ）の値が、ｌｏｇ_１０（Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）の信頼区間の外であるサンプル数の比を求めることができ、この比が一定値より大きい場合（例えば、一定値は、０．４以上の数であり得る）に、２つのＲｅｐＳｅｑが一つの細胞に由来すると推定することができる。

【0084】

また、同じ第一のクラスターに分類された２つのＲｅｐＳｅｑは同一の細胞内に存在することが分かり、したがって、この２つのＲｅｐＳｅｑは、ｃＯＴＵに分類されうる。このようにして、第二のクラスターをｃＯＴＵとして再分類することができる。

【0085】

あるいは、所定の遺伝子が、１６ｓ ｒＲＮＡである場合には、ＲＤＰ分類を用いて分類することや、ＲＤＰ分類における１６ｓ ｒＲＮＡのトレーニングセットを機械学習させることによって、最も高いスコアを有する予測分類群を作成し、これをｃＯＴＵとすることができる。なお、ＲＤＰ分類は、リボソームデータベースプロジェクトにより開発された１６Ｓ ｒＲＮＡの塩基配列から微生物種を判別するためのツールである。

【0086】

また、本発明の方法は、本発明の方法で算出された細胞の総数を、光学顕微鏡等によるカウントから推定される細胞の総数によって補正（または標準化）することをさらに含んでいてもよい。本発明の方法で算出された細胞の総数を補正（または標準化）することにより、細胞数（例えば、特定のクラスターの細胞数、または特定のｃＯＴＵの細胞数）の予測精度を向上させることができる。

【0087】

本発明の方法は、２つの異なる細胞集団間の比較を行うことに用いることができる。そして、本発明の方法は、
（Ｅ）第一の細胞集団と、第一の細胞集団とは異なる第二の細胞集団のそれぞれに関して、細胞集団に含まれる（ｉ）ｃＯＴＵの数および／または（ｉｉ）特定のｃＯＴＵに含まれる細胞の数を推定し、第一の細胞集団に関して推定された（ｉ）ｃＯＴＵの数および／または（ｉｉ）特定のｃＯＴＵに含まれる細胞の数を、第二の細胞集団に関して推定された（ｉ）ｃＯＴＵの数および／または（ｉｉ）特定のｃＯＴＵ含まれる細胞の数と比較すること
をさらに含んでいてもよい。

【0088】

上記（Ｅ）において、比較する細胞集団の細胞数を同等としておくことができる。上記（Ｅ）において、２つの異なる細胞集団の間で、ｃＯＴＵの数や各ｃＯＴＵに含まれる細胞の数を比較することによって、各細胞集団の特徴をｃＯＴＵの観点で記述することもできる。

【0089】

２つの細胞集団は、例えば、同一対象の同一部位から異なる時間に単離された細胞集団であり得、同一対象の異なる部位から同一時間に単離された細胞集団であり得、異なる対象の同一部位から同一時間に単離された細胞集団であり得る。

【0090】

同一対象の同一部位から異なる時間に単離された細胞集団を上記（Ｅ）により比較した場合には、サンプルの取得時間の相違（例えば、健康状態の変遷、治療前後の健康状態の相違、疾患や状態の発症や発展）が、ｃＯＴＵの観点で記述されることとなる。また、同一対象の異なる部位から同一時間に単離された細胞集団を上記（Ｅ）により比較した場合には、取得部位の相違（例えば、臓器毎の細菌叢の相違）が、ｃＯＴＵの観点で記述されることとなる。また、異なる対象の同一部位から同一時間に単離された細胞集団を上記（Ｅ）により比較した場合には、対象の相違（例えば、健康状態、性別、地域、人種等）が、ｃＯＴＵの観点で記述されることとなる。

【0091】

本発明の方法は、
（Ｆ）第一の細胞集団に関して推定された（ｉ）ｃＯＴＵの数および（ｉｉ’）特定のｃＯＴＵに含まれる細胞の数と、第二の細胞集団に関して推定された（ｉ）ｃＯＴＵの数および（ｉｉ’）特定のｃＯＴＵに含まれる細胞の数とを比較すること
をさらに含んでいてもよい。

【0092】

上記（Ｆ）では、第一の細胞集団について推定されたｃＯＴＵの数と第二の集団について推定されたｃＯＴＵの数との相関を決定することができる。

【0093】

上記（Ｆ）ではまた、第一の細胞集団について推定された１以上の特定のｃＯＴＵと、当該１以上のｃＯＴＵに対応する第二の細胞集団について推定された１以上のｃＯＴＵとを比較することができる。ここで、ある細胞集団から推定されたｃＯＴＵと、他の細胞集団から推定されたｃＯＴＵが対応するかどうかは、当該ｃＯＴＵに含まれるすべての塩基配列（またはエラー修正後の塩基配列）が同一であるか否かによって確認することができる。上記（Ｆ）では特に、各ｃＯＴＵに含まれる細胞数の増減が、他のｃＯＴＵに含まれる細胞数の増減と正に相関するか、負に相関するか、または相関しない（相関が弱い）かを決定することができる。これにより、ｃＯＴＵ間のネットワークを推定することができる。

【0094】

または、群集生態学の分野において群間の類似度の指標となる様々な指標を用いて細胞集団（この細胞集団は、複数のｃＯＴＵ分類群を含み、各ｃＯＴＵ分類群に対して細胞数が決定されている）を比較することができる。例えば、第一の細胞集団と第二の細胞集団との類似度は、各ｃＯＴＵに含まれる細胞数の差の２乗平均平方根（ｃ．ｆ． ユークリッド距離）として求められ得る。また、第一の細胞集団と第二の細胞集団との類似度は、各ｃＯＴＵに含まれる細胞数の差の絶対値の和（ｃ．ｆ． マンハッタン距離）として求められ得る。これらの数値は、数字が大きければ大きいほど非類似であることを示し、完全に同一である場合、０である。Ｂｒａｙ－Ｃｕｒｔｉｓ非類似度（指数）は、マンハッタン距離を標準化したものである。第一の細胞集団の細胞組成が（Ｘ_１１，・・・・，Ｘ_１ｎ）とし、第二の細胞集団の細胞組成が（Ｘ_２１，・・・・，Ｘ_２ｎ）であるとき、Ｂｒａｙ－Ｃｕｒｔｉｓ指数は以下の式により求められる。

【0095】

【数1】

【0096】

Ｂｒａｙ－Ｃｕｒｔｉｓ指数は、２つの群が全く異なる場合に１となり、完全に一致するときに０となる。このように、異なる場合に大きくなるように設計された指数であるため、非類似度と呼ばれることがある。Ｂｒａｙ－Ｃｕｒｔｉｓ指数は、統計処理プログラム（例えば、Ｒのパッケージｖｅｇａｎの関数、例えば、ｖｅｄｉｓｔ関数）を用いて計算することができる。その他、類似度は、Ｍｏｒｉｓｈｉｔａ指数、Ｊａｃｃａｒｄ指数、Ｃｈａｏ指数などの、群衆生物学の分野においてよく用いられる評価指標によって評価することができる。推定された類似度の標準偏差および信頼区間については、ブートストラップ法などによって評価できる。

【0097】

本発明の方法は、
（Ｇ）ｃＯＴＵ間の階層的クラスタリングを実施することをさらに含み得る。
階層クラスタリングは、例えば、ｃＯＴＵ間の相関の強度（例えば、スピアマンの相関係数ｒ）に基づいて、当業者に周知の方法により行うことができる。階層クラスタリングは、また、ｒから算出されるｃＯＴＵ間の距離に基づいて、当業者に周知の方法により行ってもよい。距離は例えば、１－最小（│ｒ’│）［ｒ’∈（ｒ － ＯＣＩ， ｒ＋ＯＣＩ）］｛ここで、ＯＣＩは各ｒの９０％片側信頼区間を意味する｝により算出され得る。階層的クラスタリングの結果は、系統樹として表示することができる。これは、例えば、Ｒのパッケージｈｃｌｕｓｔまたはｐｈｅａｔｍａｐにより行うことができる。また、ピアソンの相関係数ｒが閾値（例えば、０．５以上、または０．６以上など）以上となるｃＯＴＵのネットワークをパッケージｉｇｒａｐｈを用いて図示することができる。このようにして、複数の細胞集団におけるｃＯＴＵの関係性から、ｃＯＴＵ間の相関を図示することができる。

【0098】

ｃＯＴＵが既知の微生物に対応する場合には、既知の微生物間の相関を明らかにすることができるが、ｃＯＴＵが未知の微生物に対応するものであっても、ｃＯＴＵ間の相関を明らかにすることができる。あるｃＯＴＵが、既知の微生物の一つに対応する場合、既知の微生物の他の微生物（当該他の微生物は、未知であっても既知であってもよい）との相関を明らかにすることができる。また、相関する２つのｃＯＴＵが、既知の２つの微生物に対応する場合、２つの既知の微生物間の新しい相関を見出すことなどに用いることができる。このようにして、相関するｎ個のｃＯＴＵが、既知のｎ種類の微生物に対応する場合には、ｎ種類の既知の微生物間の新しい相関を見出すことができる。このように本発明の方法は、複数の細胞集団（例えば、複数の細菌叢）を調べることによって、微生物間の相関を明らかにすることに用いることができる。対象の健康状態は、その対象が有する細菌叢と相関する可能性がある。したがって、対象の健康状態とあるｃＯＴＵとの相関をさらに調べることによって、ｃＯＴＵ自体が未知の微生物であったとしても、当該未知の微生物に対応するｃＯＴＵから対象の健康状態を予測することが可能となる｛ｃＯＴＵ自体は、異なるサンプル間でも共通することに留意されたい｝。また、対象の健康状態とあるｃＯＴＵとの相関に加えて、当該ｃＯＴＵと相関する他のｃＯＴＵとの相関をさらに調べることによって、ｃＯＴＵからの対象の健康状態の予測精度を向上させることができると期待できる。

【0099】

このように、これまで１つの微生物が１つの所定の遺伝子しか有しないことを前提とした分析をしていたのに対して、本発明では、１つの微生物が複数の所定の遺伝子を有する場合であっても、それをｃＯＴＵという新しい群の概念により記述する方法を提供する。さらに、ｃＯＴＵ毎に、細胞バーコードを用いて細胞数を定性的にカウントすることにより、各ｃＯＴＵに含まれる細胞の正確な計数が可能である。分析対象となる細胞集団を本発明の方法によって分析すれば、内在するｃＯＴＵの種類と、各ｃＯＴＵに含まれる細胞数を決定することができる。取得される細胞集団と、ｃＯＴＵの種類および各ｃＯＴＵに含まれる細胞数とを分析することによって、未知の微生物を含む細胞集団に関しても、当該未知の微生物の情報を用いたより詳細な細胞集団の分析が可能となる。

【0100】

さらには、本発明は、遺伝子のコピー数が細胞によって異なる場合でも、分析精度が落ちない利点を有する。すなわち、同一種の微生物であっても、細胞内での遺伝子コピー数が細胞により異なる場合がある。このような場合、従来の方法では、遺伝子のコピー数が細胞カウントに影響する可能性があった。本発明の方法では、細胞バーコードを用いて細胞数を定性的にカウントするので、細胞内の遺伝子のコピー数には影響を受けずに細胞の計数ができ得る。微生物によっては、環境に影響する物質（例えば、毒物や増殖因子など）を放出する。細胞の数を正確に測ることは、放出される物質の量をより正確に推定することが可能となり、当該放出される物質量に基づく数理モデリングの途を拓き得る。

【0101】

本発明の方法では、シークエンスされる遺伝子は、特定の１種の遺伝子であってもよいし、複数の遺伝子であってもよい。本発明の方法では、シークエンスされる遺伝子は、全ゲノムである必要はない。

【0102】

また、従来法では、例えば、１６Ｓ ｒＲＮＡの分析において、細胞集団に含まれる全１６Ｓ ｒＲＮＡをコードする遺伝子の塩基配列をシークエンスして、得られた塩基配列を閾値に基づいて分類していた。閾値としては、例えば、同一性において９７％と設定し、９７％以上の同一性を有するものを同じ遺伝子とみなして分析した。しかし、このような分析では、異なる種、異なる属、異なる科等の本来的に生物学的に異なる分類群に属するべき微生物が１つの群として認識されることとなっていた。しかし、本発明の方法では、ある新規な１６Ｓ ｒＲＮＡが、本当に新規なものであるか、実験的なエラーによるものであるかを判別することができる。例えば、複数の細胞において発見される同一配列は、本来的に存在していた配列である可能性があるが、これを細胞バーコードによって確認することができる。このようにして、本発明の方法では、塩基配列が相違すれば、塩基配列の類似度に影響を受けない評価方法となり得る。

【実施例】

【0103】

実施例においては、既知の濃度の既知の細菌で構成した模擬的な細菌叢（ここでは「模擬細胞集団」という）を作製して測定系を検証し、その後、実際の細菌叢（ここでは盲腸の細菌叢）を調べた。

【0104】

［方法］
模擬細胞集団の準備
ヒト腸内細菌株（ＡＴＣＣ２９０９８、ＡＴＣＣ７００９２６、ＤＳＭ１４４６９、ＪＣＭ１２９７、ＪＣＭ５８２４、ＪＣＭ５８２７、ＪＣＭ９４９８、ＪＣＭ１０１８８、ＪＣＭ１４６５６、およびＪＣＭ１７４６３）からなる模擬細胞集団を調製した^９。これらの株の名称、供給源、培地および培養条件を表１に示す。培養菌を１０％グリセロールで元の培地に保存するか、リン酸緩衝生理食塩水（ＰＢＳ）中で、実験まで－８０℃で保存した（表１）。ＪＣＭ１４６５６およびＤＳＭ１４４６９は、培養後に遠心分離を用いてＰＢＳによって１回洗浄された。ＪＣＭ１０１８８をＧＡＭ寒天（ニッスイ）上で培養し、細菌コロニーを収集し、３，２００ｒｐｍで１分間ボルテックスすることによりＰＢＳに懸濁した（ＶＯＲＴＥＸ ＧＥＮＥ ２、Ｓｃｉｅｎｔｉｆｉｃ Ｉｎｄｕｓｔｒｉｅｓ）。
１０株をＰＢＳで希釈し、クラスＩＩのバイオセーフティキャビネット内で設定された濃度に従って混合した（表１）。希釈または混合の各工程に続いて、３，２００ｒｐｍで１分間ボルテックスした。この混合１０株の「模擬細胞集団」と呼ぶ。模擬細胞集団は実験まで－８０℃で保存した。

【0105】

【表1】

【0106】

上記表中記号は以下の通りである。
＊顕微鏡画像により測定された模擬細胞集団形成における添加濃度 （ｍｅａｎ ± ｓ．ｄ．， ｎ＝５， ｃｅｌｌｓ／μｌ）。
＊＊ “＋”： グラム陽性； “－”： グラム陰性．
＃ －８０℃で１０％グリセロール下で培地中に保存された．
＃＃ －８０℃でリン酸緩衝生理食塩水（ＰＢＳ）中に保存された．
ＧＡＭ、Ｇｉｆｕ Ａｎａｅｒｏｂｉｃ Ｍｅｄｉｕｍ （ニッスイ）．
ＧＡＭ Ａｇａｒ、Ｍｏｄｉｆｉｅｄ ＧＡＭ Ａｇａｒ （ニッスイ）
ＬＢ、Ｌｕｒｉａ－Ｂｅｒｔａｎｉ （Ｎａｃａｌａｉ Ｔｅｓｑｕｅ）．
ＰＹＧ、Ｐｅｐｔｏｎｅ Ｙｅａｓｔ Ｇｌｕｃｏｓｅ， ＤＳＭＺ ｍｅｄｉｕｍ １０４．
ＡＴＣＣ ｍｅｄｉｕｍ １２４９、Ｍｏｄｉｆｉｅｄ Ｂａａｒ’ｓ ｍｅｄｉｕｍ ｆｏｒ ｓｕｌｆａｔｅ ｒｅｄｕｃｅｒｓ．

【0107】

顕微鏡画像による細菌濃度測定
顕微鏡下蛍光イメージングにより各株の濃度を測定した。蛍光染色した細菌を、ポリスチレンミクロスフェア（Ｂａｃｔｅｒｉａ Ｃｏｕｎｔｉｎｇ Ｋｉｔ、Ｔｈｅｒｍｏ Ｆｉｓｈｅｒ Ｓｃｉｅｎｔｉｆｉｃ）を用いて測定した。ヨウ化プロピジウム（Ｔｈｅｒｍｏ Ｆｉｓｈｅｒ Ｓｃｉｅｎｔｉｆｉｃ）を用い、７０℃で５分間加熱して細菌を染色した。体積は、細菌計数チャンバー（ＳＬＧＣ）を用いて測定したマイクロスフェアの濃度に基づいて計算した。各菌株について、５つの独立した測定を実施した；これらの５つの測定の平均濃度および標準偏差（誤差バーとして）を、模擬細胞集団における各菌株の濃度を計算するために使用した。

【0108】

【表2】

ａ，ｂ，ｃは、３回のサンプリングの反復による。
絶対濃度は、液滴デジタルＰＣＲ（ｄｄＰＣＲ）によって測定された模擬細胞集団の総濃度（９４，４００細胞／μｌ）を用いて、シークエンスから決まった生のカウント値を正規化して得られた。

【0109】

１６Ｓ ｒＲＮＡ遺伝子のサンガーシークエンス
各株の１６Ｓ ｒＲＮＡ遺伝子を、２×ＫＡＰＡ ＨｉＦｉ Ｈｏｔ ｓｔａｒｔ ｒｅａｄｙ ｍｉｘ（Ｒｏｃｈｅ）およびプライマーＦ１－ｆｕｌｌ－Ｆｗ／Ｆ３－ｆｕｌｌ－Ｒｖ（表３）を用いて増幅した。次に、増幅した１６Ｓ ｒＲＮＡ遺伝子をｐＣＲ－Ｂｌｕｎｔ ＩＩ－ＴＯＰＯベクターにクローン化し、Ｚｅｒｏ Ｂｌｕｎｔ ＴＯＰＯ ＰＣＲ Ｃｌｏｎｉｎｇ Ｋｉｔ （Ｔｈｅｒｍｏ Ｆｉｓｈｅｒ Ｓｃｉｅｎｔｉｆｉｃ）を用いてＥ． ｃｏｌｉで増幅した。次に、Ｔ７－プロモーターおよびＳＰ６－プロモーターをプライマーとして、Ｅ． ｃｏｌｉの単一コロニーからそれぞれ１６Ｓ ｒＲＮＡ遺伝子を増幅した（表３）。最後に、各コロニーから増幅した１６Ｓ ｒＲＮＡ遺伝子のＶ３－Ｖ４領域を、Ｆ２－Ｒｖプライマー（表３）を用いてサンガーシークエンス（ＦＡＳＭＡＣ）によりシークエンスした。

【0110】

【表3】

【0111】

１６Ｓ ｒＲＮＡのシークエンス法
簡単に述べると、模擬細胞集団の細菌をＰＢＳに懸濁し、細胞溶解のためにリゾチーム、アクロモペプチダーゼ、およびプロテイナーゼＫに連続的に供した。次いで、フェノール－クロロホルム抽出によりＤＮＡを回収した。Ｉｌｌｕｍｉｎａ ａｄａｐｔｅｒ ｏｖｅｒｈａｎｇ ｎｕｃｌｅｏｔｉｄｅ ｓｅｑｕｅｎｃｅ（表３のＣＯＮＶ３４１ＦおよびＣＯＮＶ８０５Ｒ）を含む領域特異的プライマーを用いて、１６Ｓ ｒＲＮＡ遺伝子のＶ３－Ｖ４領域を増幅した。増幅産物を、ＡＭＰｕｒｅ ＸＰ磁気ビーズ（Ｂｅｃｋｍａｎ Ｃｏｕｌｔｅｒ）を用いて精製し、Ｎｅｘｔｅｒａ ＸＴ Ｉｎｄｅｘ Ｋｉｔ ｖ２（Ｉｌｌｕｍｉｎａ）を用いてインデックス化した。ＡＭＰｕｒｅ ＸＰを用いた精製後、プールしたライブラリーをＴａｐｅＳｔａｔｉｏｎ（Ａｇｉｌｅｎｔ）およびＫＡＰＡ Ｌｉｂｒａｒｙ Ｑｕａｎｔｉｆｉｃａｔｉｏｎ Ｋｉｔ ｆｏｒ Ｉｌｌｕｍｉｎａ （Ｋａｐａ Ｂｉｏｓｙｓｔｅｍｓ）により定性および定量した。２０％ＰｈｉＸ ｃｏｎｔｒｏｌ ｖ３（Ｉｌｌｕｍｉｎａ）をスパイクした変性ライブラリーを、ＭｉＳｅｑプラットフォーム（Ｉｌｌｕｍｉｎａ、２×３００ｂｐ ｐａｉｒｅｄ－ｅｎｄ ｒｅａｄｓ）でシークエンスした。配列データを質について確認し、Ｔｒｉｍｍｏｍａｔｉｃ ｖｅｒｓｉｏｎ ０．３８^４７を用いてトリミングした。ＯＴＵはＭｏｔｈｕｒバージョン１．３５．１^４８を用いて９７％の同一性閾値でクラスター化した。各ＯＴＵで最も豊富に存在する配列は、ＯＴＵの代表的な配列として選択された（図１ｂ）。

【0112】

マウスの準備
マウスの処置はすべて、理研の施設内動物実験委員会が承認したプロトコールに基づき、研究所の倫理指針に準拠して実施した。状態を維持しているマウスは以下のとおりであった。６週齢のＣ５７ＢＬ６／Ｊ雄マウスをＣＬＥＡ Ｊａｐａｎから購入し、サンプリング前に同じケージにＣＥ－２飼料（ＣＬＥＡ Ｊａｐａｎ）を給餌することにより理研施設で３日間維持した。

【0113】

マウス盲腸内容物の採取
ネズミ盲腸は、セボフルラン麻酔下で頚椎脱臼後１０分以内に手術により外に出した。異なる部位の盲腸内容物（図２ａ）を、滅菌済みはさみを用いてスライスすることによりサンプリングした。サンプリングプロセスは、クラスＩＩのバイオセーフティキャビネット内で手術後１０分以内に行われた。各マウスの各部位のサンプルを、ＤＮＡ Ｌｏｂｉｎｄ Ｔｕｂｅ（Ｅｐｐｅｎｄｏｒｆ）に採取した。対照では、２本の空試験管を用いた。サンプルの重量は、ＤＮＡ Ｌｏｂｉｎｄ Ｔｕｂｅ（全サンプルについて８．５７～１９．８２ｍｇの範囲）に採取した直後に測定した。次に、各試料を添加したＰＢＳ（１ｍｇ当たり５０μｌ）に分散させ、３，２００ｒｐｍで１分間ボルテックスすることにより混合した。懸濁試料は、その後の実験まで４℃で保存した。

【0114】

細胞外ＤＮＡの濾過
マウス盲腸試料を、１ｍｇの盲腸内容物当たり１ｍｌのＰＢＳに希釈し、その後、３，２００ｒｐｍで１分間ボルテックスした。対照では、空のチューブにＰＢＳを添加した。次いで、０．２２μｍの孔径のＵｌｔｒａｆｒｅｅ－ＭＣ Ｃｅｎｔｒｉｆｕｇａｌ Ｆｉｌｔｅｒ（Ｍｅｒｃｋ）を用いて、希釈した試料４００μｌを遠心分離（１０，０００ｇ、１０分、４℃）により濾過した。メンブラン上に残った試料に４００μｌの新鮮なＰＢＳを添加し、ピペッティングにより懸濁した後、全量を新しいＤＮＡ Ｌｏｂｉｎｄ Ｔｕｂｅ中に写した。次いで、懸濁した試料を、３，２００ｒｐｍで１分間ボルテックスした。懸濁試料およびフィルター通過液に含まれる細胞外ＤＮＡは、その後の測定まで４℃で保存した。なお、０．２２μｍのフィルターによるＤＮＡ分離の適切性は、０．１μｍの孔径のフィルターを用いた場合と比較して、フィルター通過液中の細胞外ＤＮＡの量がほぼ一致し、フィルター通過液中に細胞が検出されなかったこと、フィロターろ過後のフィルター上から回収された細菌量が等しいこと、および、フィルターから回収された細菌量がデジタルＰＣＲによる細菌数と相関することから確認された（図３１参照）。

【0115】

ＢａｒＢＩＱ法
総濃度測定
細胞または細胞外１６Ｓ ｒＲＮＡ遺伝子の総濃度は、プライマーＦ１－ＦｗおよびＦ１－Ｒｖ（表３）を用いて、Ｄｒｏｐｌｅｔ Ｄｉｇｉｔａｌ^ＴＭ ＰＣＲ（ｄｄＰＣＲ）（Ｂｉｏ－Ｒａｄ）により測定した。等モル混合した４つの細胞バーコードテンプレート（表３； ２４個のランダム塩基を含有する各テンプレートは、我々の以前の文献^２５に従って設計され、ランダム塩基の数は、単一のＭｉＳｅｑシークエンス作業において測定された個々の細胞を区別するのに十分であった）の濃度も、プライマーＮｏＢｉｏｔｉｎ－Ｌｉｎｋ－ｂａｒｃｏｄｅ－ＦおよびＰ５－ｉｎｄｅｘ－Ｒ１Ｐ－ｂａｒｃｏｄｅ－Ｒ（表３）を用いたｄｄＰＣＲによって測定した。ｄｄＰＣＲは、ＱＸ２００^ＴＭ ｄｄＰＣＲ^ＴＭ ＥｖａＧｒｅｅｎ（商標） Ｓｕｐｅｒｍｉｘ（Ｂｉｏ－Ｒａｄ）のユーザーマニュアルに従って実施した。

【0116】

ワンステップの液滴増幅
シークエンスライブラリーを作製するために、合計約２４０，０００個の細胞（または細胞外の１６Ｓ ｒＲＮＡ遺伝子の２０，０００コピー）を、等モル混合細胞バーコード、プライマー（４００ ｎＭ Ｐ７－Ｒ２Ｐ－３４１Ｆ、４００ ｎＭ Ｐ５－ｉｎｄｅｘ－Ｒ１Ｐ－Ｒ、１０ ｎＭ Ｂｉｏｔｉｎ－ｌｉｎｋ－８０５Ｒ、および１０ ｎＭ Ｂｉｏｔｉｎ－Ｌｉｎｋ－Ｆ）、ｄｄＰＣＲ^ＴＭ Ｓｕｐｅｒｍｉｘ ｆｏｒ Ｐｒｏｂｅｓ （Ｎｏ ｄＵＴＰ）（Ｂｉｏ－Ｒａｄ）、１２８ユニットのＰｌａｔｉｎｕｍ Ｔａｑ （Ｉｎｖｉｔｒｏｇｅｎ）、および１００ ｎＭ ＮＴＰを含む溶液９６０μＬと混合した。３，２００ｒｐｍで１分間ボルテックス後、混合溶液をＢｉｏ－Ｒａｄ ｄｒｏｐｌｅｔ ｇｅｎｅｒａｔｏｒにより液滴に封入し、３０μｌの混合溶液および８０μｌのＤｒｏｐｌｅｔ Ｇｅｎｅｒａｔｉｏｎ Ｏｉｌ ｆｏｒ Ｐｒｏｂｅ （Ｂｉｏ－Ｒａｄ）をＤＧ８^ＴＭカートリッジ上の各チャネルにロードした（各試料について３２チャネルを使用した）。模擬細胞集団測定のために、約６００，０００個の細胞を、細胞バーコードの約６００，０００コピー、３２０ユニットのＰｌａｔｉｍｕｍ Ｔａｑ、およびプライマー、ｄＮＴＰ、ｄｄＰＣＲ^ＴＭ Ｓｕｐｅｒｍｉｘ ｆｏｒ Ｐｒｏｂｅｓ （Ｎｏ ｄＵＴＰ）を含む２４００μＬの溶液と混合し；次いで、ボルテックス後、混合溶液を、試料あたり８０チャンネルを使用して液滴に封入した。Ｍｉｓｅｑシークエンスのためのライブラリーを、液滴中におけるワンステップＰＣＲによって生成した（９５℃の５分間；９４℃の４５秒間および６０℃の１５０秒間の６サイクル；９４℃の２５秒間および６０℃の８０秒間の４９サイクル；９８℃の１０分間）。

【0117】

ライブラリーの回収および精製
液滴増幅技術により生成したライブラリーをクロロホルムを用いて回収し、８０μｌのＴＥバッファー（Ｉｎｖｉｔｒｏｇｅｎ）および２８０μｌのクロロホルム（Ｓｉｇｍａ）を、各ＤＧ８^ＴＭカートリッジ（８ウェル）から収集した液滴と混合し、その後１０回ピペッティングし、水および有機相が分離されるまでボルテックスし；遠心後（２１，９００ｇ、１０分）、ライブラリーを含む水相の溶液を抽出した。次いで、ＡＭＰｕｒｅ ＸＰを用いたビーズ精製および２％Ｅ－Ｇｅｌ^ＴＭ ＥＸ Ａｇａｒｏｓｅ Ｇｅｌｓ （ＴｈｅｒｍｏＦｉｓｈｅｒ Ｓｃｉｅｎｃｔｉｆｉｃ）を用いたゲル精製により、未連結バーコード増幅物、残存プライマー、回収溶液中の副産物などの非標的ＤＮＡを除去した。その後、ビオチン化された結合していない１６Ｓ ｒＲＮＡ増幅物をストレプトアビジン磁気ビーズ（ＮＥＢ）により除去し、結合していない１６Ｓ ｒＲＮＡ増幅物をプライマーのＢｉｏｔｉｎ－ｌｉｎｋ－８０５Ｒ（図５）によりビオチン化した^２８。ＡＭＰｕｒｅ ＸＰ、ゲル、およびストレプトアビジンビーズを用いた精製工程は、それぞれ２回実施した。最後に、精製したライブラリーをＤＮＡ Ｃｌｅａｎ ａｎｄ Ｃｏｎｃｅｎｔｒａｔｏｒ Ｋｉｔ（Ｚｙｍｏ Ｒｅｓｅａｒｃｈ）により濃縮した。ライブラリーの品質をＡｇｉｌｅｎｔ ２１００ Ｂｉｏａｎａｌｙｚｅｒにより確認し、プライマーＰ１＿ｑＰＣＲ＿ＦｗおよびＰ２＿ｑＰＣＲ＿Ｒｖ（表３）を用いてｑＰＣＲ（ＫＡＰＡ ＳＹＢＲ Ｆａｓｔ ｑＰＣＲ ｋｉｔ， ＫＡＰＡ Ｂｉｏｓｙｓｔｅｍｓ）により濃度を測定した。ＡＭＰｕｒｅ ＸＰ、ゲル、およびストレプトアビジンビーズを使用する精製工程の詳細なプロトコールを、各製品のユーザー指示に従って実施した。

【0118】

ＭｉＳｅｑシークエンシング
サンプルのライブラリーを、ＭｉＳｅｑ ｐｌａｔｆｏｒｍ （ＭｉＳｅｑ Ｒｅａｇｅｎｔ Ｋｉｔ ｖ３、６００サイクル、Ｉｌｌｕｍｉｎａ）上で、Ｒｅａｄ １については３０サイクル、Ｉｎｄｅｘ １については２９５サイクル、Ｉｎｄｅｘ ２については８サイクル、Ｒｅａｄ ２については２９５サイクルを割り当ててペアエンドのシークエンスをした（図５）。Ｉｌｌｕｍｉｎａ Ｉｎｄｅｘ １シークエンシングプライマーを、インデックスの代わりに１６Ｓ ｒＲＮＡ配列を読み取るためにＩ１＿ｐｒｉｍｅｒ（表３）と命名されたカスタムプライマーに置き換えた。シークエンスのための配列の不均一性を維持するために、別途作製したスパイクインコントロールを試料と共にシークエンスした（図１８および１９）。より具体的には、細菌、細胞外のＤＮＡ、または細胞バーコードの総濃度は、ＱＸ２００^ＴＭ ｄｄＰＣＲ ＥｖａＧｒｅｅｎ^ＴＭ Ｓｕｐｅｒｍｉｘ（Ｂｉｏ－Ｒａｄ）の指示に従い、Ｄｒｏｐｌｅｔ Ｄｉｇｉｔａｌ^ＴＭ ＰＣＲ（ｄｄＰＣＲ）により測定した。細菌および細胞外ＤＮＡ試料については、１６Ｓ ｒＲＮＡ遺伝子のＶ１－Ｖ２領域を標的とするプライマー、Ｆ１－ＦｗおよびＦ１－Ｒｖ、または１６Ｓ ｒＲＮＡ遺伝子のＶ３－Ｖ４領域を標的とする３４１Ｆおよび８０５Ｒを用いた（表３）。細胞バーコードには、プライマーであるＢｉｏｔｉｎ－Ｌｉｎｋ－ｂａｒｃｏｄｅ－ＦおよびＰ５－ｉｎｄｅｘ－Ｒ１Ｐ－ｂａｒｃｏｄｅ－Ｒ（ｉｎｄｅｘ ＧＴＡＣＴＧＡＣ含有）を用いた（表３）。ＱＸ２００^ＴＭ ｄｄＰＣＲ^ＴＭ ＥｖａＧｒｅｅｎ^ＴＭ Ｓｕｐｅｒｍｉｘ、１μＭプライマー、１μＭ ｄＮＴＰ、およびサンプル（多重希釈、１分間３，２００ｒｐｍのボルテックス）を３０μｌの容量で混合し、混合のためにピペットで分注した。次に、Ｄｒｏｐｌｅｔ Ｇｅｎｅｒａｔｉｏｎ Ｏｉｌ ｆｏｒ ＥｖａＧｒｅｅｎ（Ｂｉｏ－Ｒａｄ）、ＤＧ８^ＴＭカートリッジ（Ｂｉｏ－Ｒａｄ）、およびＤｒｏｐｌｅｔ Ｇｅｎｅｒａｔｏｒ （Ｂｉｏ－Ｒａｄ）を用いて、混合溶液を液滴に封入した。液滴ＰＣＲは、以下のステップにより実施した。初期変性には９５℃が５分；変性には９５℃４５秒、アニーリングと伸長には６０℃１５０秒の６サイクル；変性には９５℃２５秒とアニーリングと伸長には６０℃８０秒の３９サイクル（Ｆ１－ＦｗとＦ１－Ｒｖ）または９５℃２５秒とアニーリングと伸長には６０℃８０秒の３４サイクル（プライマー３４１Ｆ／８０５Ｒ）；シグナル安定化には４℃５分間および９０℃５分間。その後、液滴の蛍光強度をＱＸ２００ Ｄｒｏｐｌｅｔ Ｒｅａｄｅｒ （Ｂｉｏ－Ｒａｄ）により測定し、ソフトウェアＱｕａｎｔａＳｏｆｔ（Ｂｉｏ－Ｒａｄ）による強度の二峰性分布の谷である閾値に基づいて陽性および陰性液滴数を決定した（図１８ａ）。最後に、サンプルの濃度は、陽性および陰性液滴の比率およびサンプルの希釈率に基づいて算出した。

【0119】

両プライマーセット、Ｆ１－Ｆｗ／Ｆ１－Ｒｖおよび３４１Ｆ／８０５Ｒを用いて、同一試料（Ｃ５７ＢＬ６／Ｊ雄マウスから得た盲腸試料）について、両細胞と細胞外のＤＮＡの合計濃度を測定し、両者の測定濃度が一致することを確認した。以下の理由から、ＢａｒＢＩＱの細菌試料の濃度測定にはプライマーＦ１－Ｆｗ／Ｆ１－Ｒｖを用いた。

【0120】

この比較のために、３４１Ｆ／８０５Ｒ（図１８ｂ）の場合の陽性および陰性の液滴の分布の間の明らかな不明瞭な分離のため、ガウスフィッティングによって陽性および陰性の液滴の割合を決定した。Ｒ ｐａｃｋａｇｅ ｍｉｘｔｏｏｌ中の関数ｎｏｒｍａｌｍｉｘＥＭにより、４つのガウス分布を用いて強度分布のピークにフィットさせた（図１８ｃ）。簡単に言えば、２つのガウス分布によるフィットは十分であり得、一方は陽性の液滴に対して、他方は陰性の液滴に対してである。しかし、データは明らかに２つ以上のガウス分布があることを示した。従って、異なる数のガウス分布による強度分布をフィットさせた。４個以上のガウス分布（６個以下を試した）を用いたとき、陽性の液滴の割合は安定であることが分かった（図１８ｄ）。このことは４個のガウス分布が強度分布を説明するのに十分であることを示唆した。陽性の液滴の割合を計算するために、このガウス分布の平均が強度二峰性分布の見かけの谷よりも大きい場合、および陰性の液滴に対しては逆の場合、陽性の液滴としてフィットしたガウス分布を仮定した。最後に、２つのプライマーセットを用いて結果間で陽性の液滴の割合を比較し、それらが細菌細胞と細胞外ＤＮＡ試料の両方について基本的に異ならないことを見出した（図１８ｅ）。プライマーＦ１－Ｆｗ／Ｆ１－Ｒｖを用いた陽性の液滴と陰性の液滴の間の分離は、３４１Ｆ／８０５Ｒ（図１８ａ，ｂ）を用いたものよりはるかに明瞭であったため、ＢａｒＢＩＱについてはＦ１－Ｆｗ／Ｆ１－Ｒｖプライマーを選択した。

【0121】

液滴調製時の細菌濃度およびバーコード濃度の調整
液滴生成のために、濃度２５０細胞／μｌの細菌を使用した。この濃度は、１つの液滴の体積が約０．８ｎｌであるので、約２０％の液滴が細菌を含有することとなる。この条件下では、ポアソン分布に従い、細菌を含有する液滴の９０％以上は１つの細菌のみを含み、他のものは２つ以上の細菌を含むこととなる。
理論的には、ＢａｒＢＩＱは、シークエンスにより決定された各ｃＯＴＵの比例濃度を総濃度を用いて正規化することによってｃＯＴＵの絶対濃度を測定し、異なる細胞バーコード濃度は各ｃＯＴＵの比例濃度を変化させないことから、細胞バーコードの濃度はＢａｒＢＩＱにおける濃度測定に影響しない。しかしながら、より高濃度の細胞バーコードは、より多くのジャンクアンプリコンを生成し、これは１６Ｓ ｒＲＮＡ配列の同定に影響するかもしれない。一方、低濃度の細胞バーコードは、細菌の検出効率を低下させるであろう。われわれは、ＢａｒＢＩＱ測定のために１００～２５０分子／μｌの範囲の細胞バーコードを使用し、その結果、８～２０％の液滴がバーコードを含むこととなった。細胞とバーコードの両方が含まれた液滴のみがシークエンスされるため、最終的に３％～１１％がシークエンスされると予測された。

【0122】

これらの濃度による細菌細胞の検出率は３％～１１％の範囲であった。異なる試料の検出率は、同じ濃度の細胞バーコードを用いても約３倍異なり、これは細胞バーコード分子の低濃度の不安定性に起因すると思われる。シークエンスにより決定されたｃＯＴＵカウントは、細胞の異なる検出率を示す反復実験間で良好な相関を示したことから、検出率は基本的に、検出されたすべてのｃＯＴＵの比例濃度の測定には影響しないことが示唆された（図１５）。

【0123】

ＢａｒＢＩＱシークエンスのスパイクインコントロール
増幅物シークエンシング^５４においてＰｈｉｘを用いてしばしば行われるように、シークエンシングにおける不均衡な塩基型を回避するために、設計されたスパイクインコントロールをライブラリーと混合し、同時にシークエンスした。スパイクインコントロールの作成の概略を図２０に示す。最初に、１７４および１７６のランダム塩基を含む２つの一本鎖ＤＮＡ（ｓｓＤＮＡ）ＳｔｄＴａｒｇｅｔ１およびＳｔｄＴａｒｇｅｔ２を、４００ｎＭの濃度で一晩、Ｔ４ ＲＮＡリガーゼ（ＮＥＢ）により連結し、次いで６５℃下で１５分間、酵素の変性工程を行った。次に、ＳｔｄＴａｒｇｅｔ１および２の連結産物から別個の設計されたランダムバーコード（ＲａｎｄｏｍＢａｒ＿ｓｔｄ１、ＲａｎｄｏｍＢａｒ＿ｓｔｄ２、ＲａｎｄｏｍＢａｒ＿ｓｔｄ３およびＲａｎｄｏｍＢａｒ＿ｓｔｄ４；図１６および表３）を含む４種類の異なるプライマーを用いて、伸長によって４種類のランダムバーコードテンプレットを作成し、１５分間の９０℃から室温へのアニーリング工程の後、Ｋｌｅｎｏｗポリメラーゼ（ＮＥＢ）を用いて伸長を行った。カラム精製後、４つの異なるインデックス化プライマー（ＲａｎｄｏｍＢａｒ＿ｓｔｄ２についてはＩｎｄｅｘ＿ＮＳＥ５０１、ＲａｎｄｏｍＢａｒ＿ｓｔｄ３についてはＩｎｄｅｘ＿ＮＳＥ５０２、ＲａｎｄｏｍＢａｒ＿ｓｔｄ４についてはＩｎｄｅｘ＿ＮＳＥ５０５、ＲａｎｄｏｍＢａｒ＿ｓｔｄ１についてはＩｎｄｅｘ＿ＮＳＥ５０６、図１６および表３）および他の末端にｃｏｍｍｏｎ ｐｒｉｍｅｒ ｓｔｄ＿Ｒ２を用い、最後の工程によって作製した伸長したテンプレートから４種類のＤＮＡ産物を増幅した。約６００塩基対を含む産物をゲル電気泳動で精製した。Ｐ１＿ｑＰＣＲ＿ＦｗおよびＰ２＿ｑＰＣＲ＿Ｒｖプライマーを用いてさらに２回ＰＣＲを行い、より多くの産物を増幅した；ＰＣＲの各ラウンドからの産物をゲル電気泳動により精製した。スパイクインコントロールは、プライマーＰ１＿ｑＰＣＲ＿ＦｗおよびＰ２＿ｑＰＣＲ＿Ｒｖを用いてｑＰＣＲにより測定した濃度に基づいて、これらの４種類の産物を等しい割合で混合することによって作製した。

【0124】

ユニークなバーコードあたりのリード数が平均６０を越えると各ｃＯＴＵの数が飽和することから、全シークエンシング実験におけるシークエンシング深度がデジタルカウントに十分であることを確認した（図２８および３４）。

【0125】

データ処理のパイプライン
Ｂａｒ配列とｃＯＴＵ（細胞型）を同定し、各ｃＯＴＵを定量するシークエンスにより得られたデータを処理するためのパイプラインを開発した。パイプラインの主な戦略は図６に示し、各ステップの詳細はＷＯ２０１８／２３５９３８Ａおよび以下に記載される通りであった。原則として、ＭｉＳｅｑからのリードは、まず細胞バーコード（Ｒｅａｄ Ｒ１）^２５を用いてクラスター化した。次に、同じ細胞バーコードに連結された１６Ｓ ｒＲＮＡ配列（Ｒｅａｄ Ｉ１およびＲ２）を、それらの配列同一性に基づいてさらにクラスター化した。各クラスター化１６Ｓ ｒＲＮＡ配列グループに対する代表的な配列（ＲｅｐＳｅｑ）を、各配列タイプに対するリード数とそれらのシークエンス品質の両方に基づいて生成した。各ＲｅｐＳｅｑのリード数とＲｅｐＳｅｑの各配列型に対するＲｅｐＳｅｑの数の両方に応じて、考えられる誤ったＲｅｐＳｅｑを複数のステップでさらに除去した（ＷＯ２０１８／２３５９３８Ａおよび図６参照）。独特のＲｅｐＳｅｑ配列型をＢａｒ配列と名付けた。次いで、Ｂａｒ配列を、同じ液滴におけるそれらの共検出頻度に基づいてｃＯＴＵにクラスター化した。もし２つ以上のＢａｒ配列が同じ液滴で頻繁に検出されたならば、それらは同じ細菌由来の複数の１６Ｓ ｒＲＮＡ遺伝子とみなし、それらを単一のｃＯＴＵにクラスター化した。次に、各ｃＯＴＵに対する細胞数を、固有の細胞バーコードの数（すなわち、バーコードクラスター）によってカウントした。各ｃＯＴＵの絶対細胞濃度は、ｄｄＰＣＲにより測定した試料の総濃度を用いて、ｃＯＴＵのシークエンスで計数された細胞を標準化することによって決定した。さらに、サンプリングおよび／または測定中にコンタミしたｃＯＴＵを対照により同定した。
パイプラインの大部分はＰｅｒｌ（バージョン５．２２．１）で書かれており、その他はソフトウェアで実施されていた。Ｒ（バージョン３．５．１）、ヌクレオチド配列クラスタライザー（バージョン０．０．７）^２５、ｂｗａ（バージョン０．７．１５）^４９。本パイプラインで使用されているＰｅｒｌのモジュールおよびＲのパッケージは表４に列挙されている。

【0126】

【表4】

【0127】

ＢａｒＢＩＱデータ処理の詳細
我々のシークエンスにおいて、Ｒ１（３０塩基）は細胞バーコードであり、Ｉ１（２９５塩基）およびＲ２（塩基）は１６Ｓ ｒＲＮＡ配列であり、Ｉ２（８塩基）は各試料をユニークに標識するインデックスであった。３回のシークエンス作業をすべて表４にまとめた。

【0128】

ステップ１：細胞バーコードに基づいたクラスター化
細胞バーコード（Ｒ１）のリードは、当初の低品質リードの欠失を除き、以前の報告（ＷＯ２０１８／２３５９３８Ａ）の通り、配列に基づいてクラスター化した。まず、広く実施されているように^４７、４つの連続した塩基からなる少なくとも１つのウインドウを含む低品質のＲ１リード（その平均スコアは１５より小さい）を除外した。シークエンスラン１、２、および３のリードの割合はそれぞれ、０．２３％、０．０５％および０．０６％であり、このプロセスによって除外された。次に、設計された細胞バーコードの最後の４つの固定塩基と一致するＲ１リードを次のステップのために選択した。サンプルとスパイクインコントロールの両方を含む同一シークエンスランに由来するすべての距離２のパラメータを有するＲ１リードを、ソフトウェア、ヌクレオチド配列－クラスタライザー^２５を用いてクラスター化した。異なるインデックスがなされたが、同じクラスターにクラスター化されたリードは除外された。得られたクラスターをＢＣｌｕｓｔｅｒと名付けた。各リードは、２つの１６Ｓ ｒＲＮＡ配列（Ｉ１とＲ２）と細胞バーコード（Ｒ１）を有した（図６）。

【0129】

ステップ２：低品質の末端とリードＩ１とＲ２のプライマー部に基づくトリミング
この段階では、すべてのリードの末端を、リードの質およびそれらのプライマー部分に基づいて、一定の位置でトリミングした。ＭｉＳｅｑシークエンスにおけるリードの塩基の品質は、一般に、リードの末端において減少し、末端においてより多くのエラーを生じさせる^５０。データ処理の次の段階では、リードの長さは同じに保つ必要があるので、われわれは均一な閾値を適用し、１回のシークエンス作業ですべてのリードの末端をトリミングした。全リードの平均品質に基づいて、シークエンスランのトリミング位置を決定した；トリミング位置を選択する規則は、連続する２つの位置の平均品質の平均が２５より低く（連続する２つの位置の平均品質の平均を使用して、シークエンスの品質の偶発的な変動を回避できる）始めたとき、リードの頭部から最初の位置を選ぶことであった。シークエンスラン１にはトリミング位置２３１（Ｉ１）と１９４（Ｒ２）、シークエンスラン２には２９４（Ｉ１）と２６７（Ｒ２）、シークエンスラン３には２７１（Ｉ１）と２３７（Ｒ２）を用いた。さらに、各リードのプライマー部分は、Ｉ１については２１塩基、Ｒ２については１７塩基である設計されたプライマーの長さに依存して直接トリミングされた。

【0130】

ステップ３：１６Ｓ ｒＲＮＡ配列（Ｉ１とＲ２）によるクラスター化
この段階では、１６Ｓ ｒＲＮＡ配列（Ｉ１およびＲ２）に基づいて各ＢＣｌｕｓｔｅｒ内のリードをクラスター化する２つのサブステップを実施した。

【0131】

ステップ３．１：置換距離によるクラスタリング
ステップ３．１において、ソフトウェアヌクレオチド配列クラスタライザーを用いて、リードＩ１およびＲ２をそれらの間の置換距離に基づいて距離３のパラメータでクラスター化し、同じＭｉＳｅｑ ＩＤを有するリードＩ１およびＲ２を、物理的に連結することによって単一のリードとして考えた。

【0132】

ステップ３．２：リードの単一位置に基づくクラスタリング
ステップ３．１は、エラーではなく真の１６Ｓ ｒＲＮＡ配列であるかもしれない非常に類似した配列を統合したため、追加のクラスター化ステップが用いられた。ステップ３．１によって生成された各サブクラスターについて、リードは、リードの特定の位置に基づいて再びクラスター化された（すべてのリードは、第一の塩基によって整列された）。この過程の論理図を図２０に示す。リード位置ごとに、塩基 （Ａ、Ｔ、Ｃ、Ｇ） の種類を含むリード数をカウントし、１番目に豊富な塩基を含むリード数に対する、２番目に豊富な塩基を含むリード数の比（Ｒａｔｉｏ２ｎｄと命名）を算出した。さらに、各リードのカウントは、この位置での塩基のシークエンスの品質スコアによって重み付けされた；規則は、スコアが１５未満の場合は０として重み付けし、一方スコアが１５以上である場合は、当該スコアを４１で割ったスコアとして重み付けした。その後、全ての位置の中で最も高いＲａｔｉｏ２ｎｄを選択し、閾値０．７５と比較した；Ｒａｔｉｏ２ｎｄ≧０．７５の場合、２番目に豊富な塩基を含むリードは、新しいサブクラスターとして元のサブクラスターから分離された。その後、両方の新しい生成サブクラスターが同じ戦略によって再度クラスター化され、全てのサブクラスターの全ての位置のＲａｔｉｏ２ｎｄが０．７５より低くなるまで繰り返しクラスタリングが実施された。最後のサブクラスターをＳＣｌｕｓｔｅｒと名付けた（図６）。液滴中の１６Ｓ ｒＲＮＡ配列（同じ細菌からの複数の１６Ｓ ｒＲＮＡ配列）の増幅効率はしばしば偏りがあった。したがって、この場合のＲａｔｉｏ２ｎｄは０．７５未満かもしれないが、両方とも真の１６Ｓ ｒＲＮＡ配列である。幸いなことに、同一の細菌由来のこれらの異なるタイプの配列（例えば、ＡおよびＢ）の増幅バイアスはランダムに起きた。例えば、時には配列Ａがより多くのリードを有し、時として配列Ｂがより多くのリードを有したので、これら２つの配列型の両方が異なる液滴から同定されたかもしれず、増幅バイアスは、細胞計数に影響を及ぼさなかった。しかしながら、両方の配列型を検出した液滴の数は、より低い閾値を使用した場合と比較して閾値０．７５を使用した場合に減少した。このことは、同じ細菌から２つの配列を同定するのに用いたステップに影響を与えるかもしれない（ステップ１２参照）。他方、０．７５より低い閾値を使用すると、誤った配列のみを含むサブクラスターを生成し得、次のプロセスに問題を生じる。従って、閾値０．７５を用いて１６Ｓ ｒＲＮＡ配列を同定したが、同じ細菌由来の両方の配列が同じ液滴で検出された場合には別の閾値０．１を用いて検出した。閾値０．１によって生成されたデータは同じ細菌由来の複数の１６Ｓ ｒＲＮＡ配列を検出するためにのみ使用され、閾値０．７５によっても同定された１６Ｓ ｒＲＮＡ配列が使用された。なお、ステップ３では置換ミスを考慮したのみであるため、挿入ミスおよび欠失ミスはすべてＳＣｌｕｓｔｅｒとしてクラスター化したが、このパイプラインの副作用は次のステップで解決した。

【0133】

ステップ４：各ＳＣクラスターの代表的な配列（ＲｅｐＳｅｑ）の作成
各ＳＣｌｕｓｔｅｒについて、リードＩ１およびＲ２の両方に対する代表的な配列（ＲｅｐＳｅｑ）を、各塩基のシークエンス品質スコアおよび各タイプの塩基の割合の両方に基づいて生成した。各タイプの塩基の比率を計算するために、各タイプの塩基についてリード数を品質スコアで重み付けした。品質スコアが１５未満の場合には０として重み付けし；スコアが１５以上の場合には、４１で除したスコアとして重み付けした。それぞれの位置について、１番目に豊富な塩基型を代表的な塩基として用いた（図６）。単一リードによるＳＣクラスターから生成されたＲｅｐＳｅｑはエラーのリスクが高いため、この段階では単一リードによるＲｅｐＳｅｑも除去した。各ＳＣｌｕｓｔｅｒのリード数は、エラーを含む誤ったＲｅｐＳｅｑと正しいＲｅｐＳｅｑを区別する重要な情報として次の段階では使用された。

【0134】

ステップ５：シフトしたＲｅｐＳｅｑの除去
この段階で、ステップ２でプライマー部分として除外したリードの頭部（Ｉ１は２１塩基、Ｒ２は１７塩基）で生じた挿入または欠失（ｉｎｄｅｌｓ）に起因するエラー型のＲｅｐＳｅｑを除去した。たとえば、ＢＣｌｕｓｔｅｒ ｘが１６Ｓ ｒＲＮＡ配列のリードを含むと仮定し、そのうちのいくつかは頭部に２個の欠失をもつ場合、プライマー部分を切り取った後に２種類のリード（ＲｅｐＳｅｑ ｉとｊ（リードには２個の欠失がある））が生じ、ＲｅｐＳｅｑ ｊがｉの左から右へ２塩基シフトしているはずである（図２１）。このエラー型をシフトしたＲｅｐＳｅｑと名付けた。
ステップ５の論理図を図２１に示すが、戦略は以下の通りである。ａ）各ＢＣｌｕｓｔｅｒで考えられるすべてのＲｅｐＳｅｑｓ型ペアを見出し、１つのＲｅｐＳｅｑ型はもう１つのＲｅｐＳｅｑ型のシフト配列でり、そのシフトが８個未満であるＲｅｐＳｅｑ型のペアのみを選択した。ｂ）シフトしたＲｅｐＳｅｑタイプ（ＡおよびＢ）の各ペアについて、より多くのＢＣｌｕｓｔｅｒで同定されたＲｅｐＳｅｑタイプを母親（ｍｏｔｈｅｒ）とし、他は可能な限りシフトと考えた。なぜなら、一般にエラーは正しいものより少ないからである。ｃ）母親とシフトの各組について、母親のリードがシフトより多いＢＣｌｕｓｔｅｒの数（Ｎｏ_{ｍｏｔｈｅｒ}）およびその反対の事例（Ｎｏ_{ｓｈｉｆｔ}）をカウントした；母親とシフトの両方を含むＢＣｌｕｓｔｅｒのみを使用した；ｄ）次に、エラーは正しいものより少ないので、Ｎｏ_{ｍｏｔｈｅｒ}がＮｏ_{ｓｈｉｆｔ}より大きい時に、母親とシフトを保存した。ｅ） 母親を伴うＢＣｌｕｓｔｅｒの中にこのシフトが存在する場合、このＢＣｌｕｓｔｅｒのシフトを削除し、母親にこのシフトのリード数（次のステップで母親に関するリード総数を使用した）を加え、母親を伴わないＢＣｌｕｓｔｅｒの中にシフトがある場合は、母親をシフトで置き換え（もし同一のシフトに関して２以上の母親が存在する場合には、より多くのＢＣｌｕｓｔｅｒで同定された母親を選ぶ）、シフトのリード数を置き換えられた母親に関するリード数として用いるという規則を用いて、ｄ）において保存された母親とシフトのペアに基づいてシフトを除去した。Ｉ１およびＲ２ＲｅｐＳｅｑは独立して処理された。

【0135】

ステップ６：Ｉ１とＲ２のＲｅｐＳｅｑの連結
このステップで、Ｉ１ ＲｅｐＳｅｑおよびＲ２ ＲｅｐＳｅｑは、それらの末端におけるそれらのオーバーラップした配列に基づいて連結された。Ｖ３－Ｖ４領域における１６Ｓ ｒＲＮＡ遺伝子の長さの分布は、Ｓｉｌｖａのデータベース（ｖ１２３．１）に依存してほぼ（＞９９．９％）４００ｂｐ～５００ｂｐの範囲であるため（図２２）、Ｉ１とＲ２の両方のリードのための２９５塩基のシークエンスは、基本的に、Ｉ１とＲ２のリードの各対の末端の間で９０以上の重複塩基を達成することができる。しかしながら、各リードの末端における低いシークエンス品質（ステップ２参照）のため、実施されたシークエンスランの最良の経験に基づいて、データ処理のために用いることができるのは、Ｉ１の約２９４塩基およびＲ２の約２６７塩基のみである。それでも、６０以上の重複塩基を検出することができる。したがって、Ｖ３－Ｖ４領域における１６Ｓ ｒＲＮＡ遺伝子の全長を得るために、Ｉ１ ＲｅｐＳｅｑとＲ２ ＲｅｐＳｅｑとの間の重複配列を見出し、それらを単一のＲｅｐＳｅｑとして連結するステップを実施した。しかし、シークエンスの質が良くなかったため、シークエンスラン１ではＩ１では２３１塩基、Ｒ２では１９４塩基のみが用いられ、したがって、オーバーラップした配列は検出されず、Ｉ１とＲ２のＲｅｐＳｅｑは連結されなかった。
一般に、Ｉ１とＲ２の両方のＲｅｐＳｅｑの末端にある数個の塩基は、偶然に同じであり得る。従って、Ｉ１ ＲｅｐＳｅｑとＲ２ ＲｅｐＳｅｑの両方の末端における５個以上の同一塩基の閾値を、偽のオーバーラップを回避するために重なりとして使用したためと考えた。理論的には、偶発的な重複の可能性は（１／４）^ｂであり、ここで、ｂは重複した塩基の数であり、５塩基の事故重複の可能性は（１／４）^５≦０．０００９８である。
さらに、模擬細胞集団およびＭ０のデータでは、すべての偶発的な重複が＜５塩基であった（短いリードが使用されたため、重複は見出せない）。
重複部分のＩ１ ＲｅｐＳｅｑとＲ２ ＲｅｐＳｅｑ間の置換の違いが、リード末端部分の品質が相対的に低かったために、稀であっても起こるかもしれない。従って、これらのエラーを除去するためにもう一つのプロセスを適用した。この戦略は、ａ）上述の連結プロセスの後、連結していないＲｅｐＳｅｑを発見した。ｂ）次に、同一のＢＣｌｕｓｔｅｒ内の他のＲｅｐＳｅｑ（それぞれＩ１およびＲ２ ＲｅｐＳｅｑを直接比較）で各ＲｅｐＳｅｑを比較し、その１つの塩基が異なるＲｅｐＳｅｑを見つけた。ｃ） １塩基相違のＲｅｐＳｅｑを連結した場合、連結されていないＲｅｐＳｅｑを削除し、連結されたＲｅｐＳｅｑにそのリードを追加した。

【0136】

ステップ７： １つの挿入と欠失（１－ｉｎｄｅｌ）ＲｅｐＳｅｑの除去
このステップでは、リードの主要部分における｛リードの頭部（すなわち、プライマー部分、ステップ５参照）においてではない｝１つの挿入または欠失（１－ｉｎｄｅｌ）エラーに起因して生じたエラータイプのＲｅｐＳｅｑを除去した。ステップ３におけるクラスタリングは上述した置換のみに基づいていたため、インデルを含む誤ったリードをすべて分離し、個々のＲｅｐＳｅｑを作成した。一般に、シークエンスのリードの中ごろに起こるインデルは非常にまれであるので（Ｓｃｈｉｒｍｅｒ Ｍ ｅｔ ａｌ．， ＢＭＣ Ｂｉｏｉｎｆｏｒｍａｔｉｃｓ ２０１６； １７：１２５）この段階では１－ｉｎｄｅｌだけでなく、ステップ９で置換を伴う２塩基インデルと１－ｉｎｄｅｌだけを考えた（後述）。
ステップ７の論理図を図２３に示す。戦略は以下の通りである。ａ）各ＢＣｌｕｓｔｅｒにおいて、１－ｉｎｄｅｌの差異を持つＲｅｐＳｅｑタイプの可能なペアをすべて見出した。ｂ） エラーは一般に正しいものより少ないので、ＲｅｐＳｅｑタイプ（ＡおよびＢ）の各１－ｉｎｄｅｌ対について、より多くのＢＣｌｕｓｔｅｒで同定されるＲｅｐＳｅｑタイプを母親とし、および他のものは１－ｉｎｄｅｌであると考えた。ｃ）母親と１－ｉｎｄｅｌの各ペアについて、母親のリードが１－ｉｎｄｅｌ（Ｎｏ_{ｍｏｔｈｅｒ}）より多いＢＣｌｕｓｔｅｒ数と、その反対の事例（Ｎｏ_{１－ｉｎｄｅｌ}）を計数した；母親と１－ｉｎｄｅｌの両方を含むＢＣｌｕｓｔｅｒのみを使用した。ｄ）われわれは、Ｎｏ_{ｍｏｔｈｅｒ}がＮｏ_{１－ｉｎｄｅｌ}より大きいときにその母親と１－ｉｎｄｅｌのペアだけを残した。なぜなら、エラーのリードは一般に正しいものより少ないからである。ｅ）可能性のある１－ｉｎｄｅｌ及びその可能性のある母親を含むＢＣｌｕｓｔｅｒの数と可能な１－ｉｎｄｅｌを含む全てのＢＣｌｕｓｔｅｒの数（Ｎｏ_{１－ｉｎｄｅｌ}）との比（Ｒｓ）を計算し、条件付き文言Ｒｓ≦（Ｎｏ_{１－ｉｎｄｅｌ}－３）／ Ｎｏ_{１－ｉｎｄｅｌ}が真である場合、可能な母親と１－ｉｎｄｅｌのペアを選択した。ｆ）選択された母親と１－ｉｎｄｅｌのペアに基づき、母親とＢＣｌｕｓｔｅｒの中に１－ｉｎｄｅｌが存在する場合は、このＢＣｌｕｓｔｅｒでは１－ｉｎｄｅｌを削除し、母親に１－ｉｎｄｅｌのリード数を加え（次のステップでは、母親としてリード総数を用いる）、もし、母親を伴わないＢＣｌｕｓｔｅｒに１－ｉｎｄｅｌが存在する場合には、１－ｉｎｄｅｌのリード数を母親に置き換え（同じ１－ｉｎｄｅｌに対して母親が２以上ある場合は、より多くのＢＣｌｕｓｔｅｒにおいて同定された母親を選択する）、１－ｉｎｄｅｌのリード数を置き換えられた母親に関するリード数として用いる。

【0137】

ステップ８：キメラの除去
この段階で、キメラ増幅により生じたエラー型のＲｅｐＳｅｑキメラを除去した。キメラは常にＰＣＲ中に起こり、産物をより複雑にする。特に１６Ｓ ｒＲＮＡ増幅物の測定ではＲｅｐＳｅｑキメラは非常によく起こる^２７。
キメラを除去する論理図は図２４に示し、その戦略は以下の通りであった。ａ）各ＢＣｌｕｓｔｅｒにおいて、ＲｅｐＳｅｑタイプ（Ａ、Ｂ、Ｃ）のすべての可能な順序のキメラをチェックした；Ａの頭部がＢの頭部部分と同じであり、Ａの他の部分（Ｂを伴う）がＣの末端部分と同一であり、かつ、Ａのリード数が３つ中最大では無い場合、Ａはキメラと考えられ、ＢとＣはこのキメラの親であると考えられた。ｂ）同定されたキメラそれぞれについて、キメラを含むＢＣｌｕｓｔｅｒの数（Ｔｏｔａｌ＿Ｎｏ）およびキメラのみを含むが親は含まないＢＣｌｕｓｔｅｒの数（Ｎｏ＿ｄ）を計数した。ｃ）条件付き文言Ｒａｔｉｏ＿ｄ（＝ Ｎ＿ｄ／Ｔｏｔａｌ＿Ｎｏ）≦０．１、かつ、Ｒａｔｉｏ＿ｄ≦１／Ｔｏｔａｌ＿Ｎｏが真のとき、ＲｅｐＳｅｑｓからキメラ候補を除外した。
ＢａｒＢＩＱはキメラを１～５％しか持たず、これは従来の方法によるもの（～７０％）^２７よりはるかに低く、この工程によりキメラを除去することができたことがわかる。ＢａｒＢＩＱでキメラがほとんど生成されなかった理由は、バーコードおよびシークエンスアダプターが、分離された空間（すなわち、液滴）におけるワンステップ増幅によって、単一の細菌由来の１６Ｓ ｒＲＮＡ遺伝子に付着したことであり、これは、異なる細菌由来の１６Ｓ ｒＲＮＡ増幅物が混合されなかったことを意味する。このアプローチは、液滴およびバーコードを用いたハイスループット１６Ｓ ｒＲＮＡ遺伝子シークエンスに関する最近の研究（Ｂｏｒｇｓｔｒｏｍ Ｅ ｅｔ ａｌ．， Ｎａｔ Ｃｏｍｍｕｎ ２０１５； ６： ７１７３およびＳｈｅｔｈ ＲＵ ｅｔ ａｌ．， Ｎａｔ Ｂｉｏｔｅｃｈｎｏｌ ２０１９； ３７（８）： ８７７－８８３）でさえ実施されていない。

【0138】

ステップ９：稀なエラーＲｅｐｓｅｑの除去
このステップでは、１インデルおよび１置換エラー（ＣａｓｅＡと命名）、１インデルおよび２置換（ＣａｓｅＢと命名）、ならびに２インデル（ＣａｓｅＣと命名）を有するＲｅｐＳｅｑのような高レベルのエラーは除去された。すでに述べたように、インデルのエラーのみが、ステップ３における我々のクラスタリング方法によって生じ得るので、ここで考察する高レベルのエラーは、インデルを含む。一方、もっと複雑なエラーはきわめてまれに起こり、ステップ１０で取り除かれる。
ステップ９の論理図を図２５に示し、戦略は以下の通りである。ａ）上記相違（Ｃａｓｅ Ａ、Ｂ、Ｃ）の何れかを有する各ＢＣｌｕｓｔｅｒのＲｅｐＳｅｑタイプの可能なペアをすべて同定し、ｂ）各同定ペアのＲｅｐＳｅｑのリード数を比較した。ＲｅｐＳｅｑ（小型／大型）間のリード数の比が閾値０．２よりも低ければ、リード数の少ないＲｅｐＳｅｑを除外し、リード数を他のペアに追加した。

【0139】

ステップ１０：低カウントＲｅｐＳｅｑの除去
大半のエラーが上記のステップで除去された後も、未知のＲｅｐＳｅｑ（Ｓａｎ配列とは異なる）は依然として模擬細胞集団のデータに残っていた。しかし、いずれも少数であった。そこで、残ったＲｅｐＳｅｑｓの種類ごとにＢＣｌｕｓｔｅｒ数をカウントした。低カウントを原因とするばらつきが大きかったため、ＲｅｐＳｅｑタイプごとにサンプリングの反復（異なるサンプリングによる同一サンプルのシークエンス）に基づく平均カウントを用いた。各反復について、各ＲｅｐＳｅｑタイプのカウントついては、すべての反復の中で最高の総カウントに対するすべてのＲｅｐＳｅｑタイプの総カウントによって正規化した。次いで、各ＲｅｐＳｅｑタイプの平均カウントを、全ての反復実験から計算した。模擬細胞集団に対して３回のサンプリングを行い、３回の反復実験から平均カウントを得た。最後に、平均カウントが２未満の場合、ＲｅｐＳｅｑタイプを除外した。
このステップの後、模擬細胞集団のデータについては、Ｓａｎ配列が一致するＲｅｐＳｅｑｓタイプを除き、残ったＲｅｐＳｅｑタイプはすべて、ＰＣＲによる１塩基エラー（ステップ１１参照）またはコンタミネーション（ステップ１４参照）として合理的に説明することができる。
１反復のみまたは２反復を用いても試験し、閾値＜６（１反復）または閾値＜３（２反復）が模擬細胞集団データに対して機能することを見出した。しかし、無作為性のため、１回と２回のサンプリングは、３回のサンプリングよりもリスクが高い可能性があるため、１回と２回のサンプリングを盲腸サンプルに用いた場合の閾値としてそれぞれ＜１０と＜５を使用することとした。

【0140】

ステップ１１：１塩基エラーＲｅｐＳｅｑの除去
この段階で、ＰＣＲによって生じたと思われるＲｅｐＳｅｑタイプの一塩基エラーを除去した。このＲｅｐＳｅｑの特徴を明らかにするため、まず、各Ｓａｎ配列と１塩基またはゼロ塩基の差異を有する残存するＲｅｐＳｅｑタイプをグループに分類した（この分析に関しては、低カウントＲｅｐＳｅｑタイプを維持した、ステップ１０を参照）。次に、各グループにおける全ＲｅｐＳｅｑタイプの平均カウントの分布（図２６ａ）をプロットし、同一グループ内のＳａｎ－配列一致型ＲｅｐＳｅｑタイプの平均カウントに対する１塩基の異なるＲｅｐＳｅｑタイプの最高平均カウントの比（最高比率）を算出した（図２６ｂ）。我々は、２つのカテゴリー（図２６ｂのカテゴリー１と２）：カテゴリー１はＳａｎ－配列一致型 ＲｅｐＳｅｑ ｔｙｐｅが１，０００カウント以上のグループ、カテゴリー２はＳａｎ－配列一致型 ＲｅｐＳｅｑ ｔｙｐｅが１，０００カウント未満のグループを見出した。カテゴリー１に関して、１塩基の異なるＲｅｐＳｅｑタイプの最高平均数は２より大きく、それらの最高比率は異なるグループ間で一貫していた。これらの１塩基の異なるＲｅｐＳｅｑ型は、おそらくＰＣＲによって生じた誤りであると結論づけた。それは、実際の１６Ｓ ｒＲＮＡ配列に対する他の実際の１６Ｓ ｒＲＮＡ配列の数の比は通常、１６Ｓ ｒＲＮＡ配列の各タイプで異なるからである。そこで、Ｓａｎ－配列一致型 ＲｅｐＳｅｑタイプのカウントに対する１塩基の異なるＲｅｐＳｅｑ型のカウントの比が１／４００未満であるという閾値を用いて、これらのＲｅｐＳｅｑ型を除去するプロセスを適用した（図２６ｂ）。１反復のみが実施された場合は、データに対して１／１００の閾値を使用した。カテゴリー２では、１塩基の異なるＲｅｐＳｅｑタイプの最高平均カウントは異なるグループ間で類似しており、＜２であったが、これはこれらの１６Ｓ ｒＲＮＡ配列の低濃度に起因する可能性があり、エラーはランダムに発生し、全てのエラー配列は一致しなかった。ＲｅｐＳｅｑｓの低カウントのリスクが高いため、ステップ１０でＲｅｐＳｅｑを除外した。模擬細胞集団データにおいて、１塩基の差のみを示す２つのＲｅｐＳｅｑタイプのカウント数の間の比が１／５０以上であれば、両方のＲｅｐＳｅｑタイプがＳａｎ配列と一致することを確認した。しかし、この比率が１／４００～１／５０の範囲にある例は見つからなかった。さらに、我々の模擬細胞集団のデータでは、１つの奇妙なＲｅｐＳｅｑタイプが検出された。この配列を確認することにより、それはサン配列ＪＣＭ５８２４－ＡおよびＪＣＭ５８２４－Ｂの中央に一致するが、ＪＣＭ５８２４－Ａ／ＢのＶ３－Ｖ４領域の全長よりもはるかに短いことを見出した。ＪＣＭ５８２４－Ａ／Ｂの中央の６ｍｅｒは、１６Ｓ ｒＲＮＡ遺伝子を増幅するために使用したフォワードプライマーの３’末端と同じであり、この奇妙な配列は同一液滴中のＪＣＭ５８２４－Ａおよび／またはＪＣＭ５８２４－Ｂの全長Ｖ３－Ｖ４領域と常に同時検出され、そのカウントは常に非常に稀であった（３回の反復において２／４／１）ことから、この奇妙なＲｅｐＳｅｑ型は、ＪＣＭ５８２４からの１６Ｓ ｒＲＮＡ遺伝子の非特異的増幅産物であると解釈した。しかし、この種の短い増幅産物は盲腸検体では見つからなかったため、我々の最終パイプラインにはこれらの短い増幅産物を検出するためのいかなるステップも含めなかった。上記の全ステップ後、残ったＲｅｐＳｅｑタイプ（固有ＲｅｐＳｅｑｓ）をＢａｒＢＩＱ同定配列（Ｂａｒ配列）と名付け、各々ＩＤ番号でラベルした。

【0141】

ステップ１２：Ｂａｒ配列のｃＯＴＵへのクラスター化
ＢａｒＢＩＱの大きな利点を生かし、この段階での細胞バーコードに基づいて、同じ細菌から複数の１６Ｓ ｒＲＮＡ配列を同定した。

【0142】

このためには、２つの可能性を考慮すべきである。１つは同一液滴中の異なる細菌の混入の可能性であり、もう１つは同一細菌細胞からの異なる配列に対する増幅バイアスにより一方の配列しか検出できない可能性である。１つ目の場合はポアソン分布に依存しており、液滴発生に低濃度の細菌を使用したため極めて稀なものである。２つ目の場合は細菌濃度に影響されないものである。実験的に、液滴数に対する細菌数の比率を２０％とすることにより、これら２つの可能性を分けることができることを見出した。

【0143】

これらの２つの可能性を区別するために、著者らはＢａｒ配列の全ての可能なペアをチェックした；各ペアについて（ＢＳ＿ＡおよびＢＳ＿Ｂとして標識した）、それらの両方を含む液滴の数（Ｏｖｅｒｌａｐとして命名）、ＢＳ＿Ａのみを含む液滴の数（Ａとして命名）、およびＢＳ＿Ｂのみを含む液滴の数（Ｂとして命名）をそれぞれ計数した。これらのカウントは、上記ステップ３．２のパラメータ０．１を使用して処理されたデータに基づいている。

【0144】

理論的には、１つのペアのＢａｒ配列が異なる細菌に由来する場合、両方のＢａｒ配列が検出される液滴の数はポアソン分布に従うはずであり、同時検出された液滴の推定数（Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐと命名）は以下のように計算できる：
Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ＝（Ａ×Ｂ×μ）／液滴総数
｛ここで、液滴総数は、細胞バーコードを含む液滴の総数であり；μは、定数であり、ＰＣＲ増幅効率、シークエンス深度効果などを含み得る液滴における検出効率のための統合パラメータである｝。他方、もしＢａｒ配列が同一細菌由来であれば、両方のＢａｒ配列が検出される液滴の数はポアソン分布に従わないであろう。

【0145】

次に、ｌｏｇ_１０変換を用いてパラメータを２つの項に分けた。
ｌｏｇ_１０ （Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）＝ ｌｏｇ_１０ （Ａ×Ｂ） － ｌｏｇ_１０ （液滴総数／μ）

【0146】

第１項のパラメータＡおよびＢはデータから得ることができるが、第２項のパラメータ液滴総数およびμは個々に測定することができない。μは異なるＢａｒ配列ペアに対しても同じであると仮定し、そして、ｌｏｇ_１０ （液滴総数／μ）は各実験の全Ｂａｒ配列ペアに対して一定であると仮定した；この用語は操作上の液滴（ＯＤ）と命名した。次に、モデルｙ＝ｘ－ＯＤを用いてｌｏｇ_１０ （Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）のｌｏｇ_１０ （Ａ×Ｂ）に対するランニング中央値を当てはめてＯＤを推定した。一般に、我々のデータでは、ほとんどのＢａｒ配列のペアは異なる細菌由来であり、それらの測定Ｏｖｅｒｌａｐは理論上のＰｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐと類似していた。したがって、ｌｏｇ_１０ （Ｏｖｅｒｌａｐ）のランニング中央値｛ここで、ランニング中央値とは、一定の大きさのウインドウａの領域における中央値と、当該領域を一定の大きさのオーバーラップｂだけずらしてさらに取得される中央値と、この操作を繰り返して得られるさらなる中央値からなる中央値の群であり、ａ＞ｂである｝を用いて、ｌｏｇ_１０ （Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）のランニング中央値を模倣した。ｌｏｇ_１０ （Ｏｖｅｒｌａｐ）のランニング中央値は、ｌｏｇ_１０ （Ａ×Ｂ）に基づいて０．４のウインドウと０．２のオーバーラップで求め、０を超える中央値のみを用いた（図２７ａの赤色の白丸）。

【0147】

フィッティングさせてＯＤを得た後、ｌｏｇ_１０ （Ａ×Ｂ）＋ＯＤに対してｌｏｇ_１０ （Ｏｖｅｒｌａｐ）でデータを再プロットした（図２７ｂ）。これは、実際には、ｌｏｇ_１０ （Ｏｖｅｒｌａｐｓ）とｌｏｇ_１０ （Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）の間の関係であった。したがって、ペアが異なる細菌由来のｌｏｇ_１０ （Ｏｖｅｒｌａｐ）のデータは、ｙ＝ｘの直線にあるはずである。しかし、ノイズのためにデータは広く分布した。

【0148】

次に、シミュレーションを行い、ｌｏｇ_１０ （Ａ×Ｂ）＋ＯＤの異なる値に対するｌｏｇ_１０ （Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）の可能な分布を推定した。最初に、Ａ、ＢおよびＯＤの異なる値について、ｌｏｇ_１０ （Ａ×Ｂ）＋ＯＤの値が同じ場合には、ｌｏｇ_１０ （Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）の分布がわずかに異なり、ＡがＢと等しい場合には分布は最も広くなり；ｌｏｇ_１０ （Ａ×Ｂ）＋ＯＤの異なる値ではｌｏｇ_１０ （Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）の分布が異なることを確認した。そこで、１～１５００の範囲で、ＡおよびＢの各可能な値（Ａ＝Ｂ、整数）に対するｌｏｇ_１０ （Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）の分布とＯＤ＝ｌｏｇ_１０（５０００）の固定値を５０００００回反復してシミュレートした。ここで、Ａ＝Ｂのときに、ポワソン分布が最も広くなると考えら、その場合においてもポワソン分布に従わない配列のペアは、異なる液滴から得られた配列である可能性がより高くなると推定できるため、ここではＡ＝Ｂでシミュレーションを行った。２つのシミュレーション数間のｌｏｇ_１０ （Ａ×Ｂ）＋ＯＤの値について、より高い近接シミュレーション値の同じ分布を用いた。その後、各分布の片側信頼区間０．９９９を算出した（図２７ｂの緑線）。

【0149】

模擬細胞集団のデータについては、同一細菌由来のＢａｒ配列ペアのｌｏｇ_１０（Ｏｖｅｒｌａｐ）のすべての値が、上側０．９９９の片側信頼区間（ＵＰ９９９）よりも大きかったが、異なる細菌由来のペアの値は片側信頼区間ＵＰ９９９と同じかそれよりも小さかった（図２７ｂ、検出効率は安定していなかったため、いくつかのＢａｒ配列のペアは下側０．９９９片側信頼区間よりも低かったが、この目標には影響しないことに留意されたい）。このデータは、同一細菌由来のペアのｌｏｇ_１０（Ｏｖｅｒｌａｐ）の値は、液滴数に対する細菌数の２０％比を用いたとき、ｌｏｇ_１０（Ｐｏｉｓｓｉｏｎ＿Ｏｖｅｒｌａｐ）よりも有意に大きく、ＵＰ９９９によって容易に区別できることを示唆した。

【0150】

次に、盲腸試料の測定と同じ方法を用いてＭ０データを分析した。ＵＰ９９９周辺では、ｌｏｇ_１０（Ｏｖｅｒｌａｐ）のｌｏｇ_１０（Ａ×Ｂ）＋ＯＤに対するプロットに明確なギャップは認められなかった（図２７ｃ）。盲腸試料の良好な閾値を見出すために、公開データベースＳｉｌｖａに基づく各Ｂａｒ配列のマッピング情報を使用した。たいていのＢａｒ配列の名前はデータベースに基づいて種レベルでは決定できず、時にはより高いレベルでも決定できないので、異なる名前にマッピングされたＢａｒ配列だけに焦点を当てた。Ｍ０データの各サンプリング反復実験において、データベース中の異なる名前にマッピングされたいくつかのＢａｒ配列ペアのｌｏｇ_１０（Ｏｖｅｒｌａｐ）は、ＵＰ９９９（図２７ｃの黒丸）より大きかった。次に、別の２回のサンプリングでこれらのＢａｒ配列のペアを調べたところ、ｌｏｇ_１０（Ｏｖｅｒｌａｐ）の値はすべてＵＰ９９９よりも低かった。１回の測定で２０，０００以上のＢａｒ配列のペアがあり、０．９９９片側信頼区間から外れていたことが合理的であったため、これらのケースはアクシデントによって統計学的に生じた可能性がある。

【0151】

統計学的にまれなケースを避けるために、複数の反復実験を用い、これら２つのＢａｒ配列が同じ細菌に由来するかどうかを決定した。理論的には、同じ細菌のＢａｒ配列の結果は、異なるサンプルで同じであるべきであり、そのため、全てのサンプルをこの目的のための反復として使用することができる。次に、マウスＭａ、ＭｂおよびＭｃの盲腸由来のすべての細胞サンプルを用い、Ｂａｒ配列ペアのｌｏｇ_１０（Ｏｖｅｒｌａｐ）がＵＰ９９９より大きいことを示したサンプル数の、両方のＢａｒ配列が検出されたサンプル総数に対する比を分析した。この比をＲａｔｉｏ＿Ｐｏｓｉｔｉｖｅとよぶ。サンプル数よりもむしろ比率を用いるのは、一部のＢａｒ配列がサンプルの一部においてのみ検出され、各Ｂａｒ配列ペアに対して用いることができるサンプルの数が異なり得るためである。信頼性を担保するために、少なくとも２つのサンプルで検出されたＢａｒ配列ペアのみを用いた。さらに、いくつかのサンプルはＯＤのフィッティングが不良であることを見いだし、フィッティングによるＯＤの標準誤差が０．０８より小さいサンプルのみを選択した。Ｂａｒ配列のマッピング名に基づき、異なる名前にマッピングしたすべてのＢａｒ配列のペアはＲａｔｉｏ＿Ｐｏｓｉｔｉｖｅが低く（図２７ｄ）、分布は指数関数的に減衰し、低い可能性でしか生じていないことが示唆された。したがって、Ｒａｔｉｏ＿Ｐｏｓｉｔｉｖｅ＞０．５の閾値を用いて、同一細菌由来のＢａｒ配列ペアを同定した。

【0152】

次に、同一細菌由来の同定されたＢａｒ配列ペアに基づき、すべてのＢａｒ配列をグループに分類した。各グループは１つのＢａｒ配列または複数のＢａｒ配列を有することができる。我々は、これらのグループを細胞ベースの操作上分類単位（ｃＯＴＵ）と名付けた。この分類の戦略は、もしこのＢａｒ配列および複数のＢａｒ配列のうちの少なくとも１つが同じ細菌のものであれば、それぞれのＢａｒ配列をグループにまとめることであった。いくつかのｃＯＴＵ内で、いくつかのＢａｒ配列ペアは、上記のプロセスによって検出されなかったが、これは、液滴が２つ以上の配列を含む場合に、検出効率が低かったためであると考えられる。

【0153】

ステップ１３：各ｃＯＴＵのカウントセル数
同一のＢＣｌｕｓｔｅｒで検出されたＲｅｐＳｅｑは、同一のｃＯＴＵに属する場合、単一細胞とみなした。次いで、各ｃＯＴＵの細胞数を、細胞バーコード（ＢＣｌｕｓｔｅｒの数）に基づいて計数した。ステップ３．２でパラメータ０．７５で処理したデータを細胞数の計数に用いた。

【0154】

ステップ１４：異物混入したｃＯＴＵの除去
この段階で、コントロールに基づいて異物混入したｃＯＴＵを除去した。異物混入したｃＯＴＵを同定するために、模擬細胞集団の対照サンプルＭ０またはマウスＭａ、ＭｂおよびＭｃの盲腸サンプルの空試験管対照を用いて、同様の時間（数日間）内に同じ条件下で測定した。

【0155】

異物混入したｃＯＴＵを検出するための戦略は、以下であった：各対照について、試料中で同定されたｃＯＴＵのＢＣｌｕｓｔｅｒ数をカウントし、対照における各ｃＯＴＵのカウント数を、試料中の同一ｃＯＴＵのカウント数と比較した。模擬細胞集団の実験では、模擬細胞集団とＭ０試料のライブラリーを調製するために異なる数の液滴を使用したため、推定総液滴数によってもカウントを正規化した。他の実験では、コントロールは空のチューブで、すべての実験で同じ数の液滴を使用したため、カウントに標準化は適用しなかった。

【0156】

模擬細胞集団に関して、３つの異なるカテゴリー（Ｉ、ＩＩ、およびＩＩＩ）（図２８）を見出した。（Ｉ）対照（すなわち、Ｍ０）におけるｃＯＴＵのカウントは試料（すなわち、模擬細胞集団）よりもはるかに大きく、Ｓａｎ－配列と一致しなかった。（ＩＩ）ｃＯＴＵの数はサンプルと対照の間で同等であり（それらの平均±ＳＤは重複した）、それらもＳａｎ配列と一致しなかった。（ＩＩＩ）試料中のｃＯＴＵの数は対照よりもはるかに高く、それらはＳａｎ配列と一致した。カテゴリーＩは、対照から試料への交差混入、または環境からの混入細菌が対照における細菌と同じであること（われわれは模擬細胞集団測定のコントロールとして実際のサンプルを使用したため）、さらに、このｃＯＴＵの配列は、試料中のＳａｎ配列と一致するＢａｒ配列（≦８６％の同一性）と非常に異なり、キメラとして説明できなかった。このことは、このｃＯＴＵがキメラから生じたエラーではないことを示唆する。カテゴリーＩＩは、サンプルと対照の両方について環境からの異物混入として説明できる。なぜなら、それらの数は異なるサンプルと対照で同様であったからである。カテゴリーＩとカテゴリーＩＩの両方について、我々はｃＯＴＵをサンプルから除去した。それは、それらがおそらく環境または他のサンプルからの異物混入であったからである。カテゴリーＩＩＩは交差異物混入の可能性もあるが、試料から対照への混入、または偶然、環境からの混入細菌は試料中の細菌と同じであった。この場合、試料中のこのｃＯＴＵのカウントから対照中のこのｃＯＴＵのカウントを差し引いた値を、試料中のこのｃＯＴＵの最終細胞数として使用した。交差異物混入の可能性は、試料中のｃＯＴＵが高濃度である場合にのみ生じたが、異物混入した試料中のそれらカウントは非常に稀であった。

【0157】

マウスＭａ、ＭｂおよびＭｃのデータについては、各試料に対して１回のみの測定を行ったため、反復ＳＤの代わりにポアソンサンプリングノイズに基づくエラーバーとしてカウントの平方根を用いた。

【0158】

マウスＭａ、ＭｂおよびＭｃのデータについて、対照として２本の空の試験管を用いた。この場合、試験管２本は反復をサンプリングするよりは実験的な繰り返しであり、ポアソン分布に従わない。また、反復回数が少ないことによる事故を回避するために、対照に対してエラーバーとして３．２７×ＳＤを使用した；さらに、３．２７×ＳＤがカウントの平均の１０％より小さい場合は、エラーバーとして平均の１０％を使用した。これらのサンプルの異物混入したｃＯＴＵを除去する規則は以下の通りであった。対照のカウント＋エラーバーがサンプルのカウント－エラーバーよりも高い場合は、サンプルからこのｃＯＴＵを除去し、コントロールのカウント＋エラーバーがサンプルのカウント－エラーバーよりも低い場合は、サンプルのカウント－コントロールのカウントをサンプルのｃＯＴＵの最終カウントとして使用した。

【0159】

異物混入したｃＯＴＵの細胞数は、模擬細胞集団の測定で検出された全細胞数の約０．５％であり、Ｍａ、Ｍｂ、およびＭｃの細胞試料測定で約４％であった。

【0160】

ステップ１５：細胞濃度の算出
各ｃＯＴＵの絶対細胞濃度は、液滴デジタルＰＣＲによって測定した総濃度を用いて、ステップ１３で得られたカウントを正規化することによって算出した。

【0161】

１６Ｓ ｒＲＮＡ遺伝子データベースとの比較
３種類の公的データベース、ＧｒｅｅｎＧｅｎｅ（リリース１３＿５）^１０、Ｒｉｂｏｓｏｍａｌ Ｄａｔａｂａｓｅ Ｐｒｏｊｅｃｔ（リリース１１．５）^１１、およびＳｉｌｖａ（リリース１３１．１）^１２において同定されたＢａｒ配列と最も近い（すなわち、最も高い同一性）１６Ｓ ｒＲＮＡ遺伝子の間の配列同一性は、ＮＣＢＩ ｂｌａｓｔ（バージョン２．７．１）^５１を用いて算出した。

【0162】

ＲＤＰ分類による分類学的予測
同定されたｃＯＴＵの門から属への分類を、ブートストラップカットオフ５０％^３６を用いたＲＤＰ分類によるそれらのＢａｒ配列に基づいて予測した。ＲＤＰ分類は１６Ｓ ｒＲＮＡトレーニングセット^１１（ｈｔｔｐｓ：／／ｒｄｐ．ｃｍｅ．ｍｓｕ．ｅｄｕ／ｃｌａｓｓｉｆｉｅｒ／ｃｌａｓｓｉｆｉｅｒ．ｊｓｐ）によりトレーニングされた。複数のＢａｒ配列を含むｃＯＴＵに対して、最も高いスコアを有する予測分類群を選択した。

【0163】

Ｂｒａｙ－Ｃｕｒｔｉｓの非類似度
Ｒパッケージｖｅｇａｎのｖｅｇｄｉｓｔ関数を用いて、細胞濃度に基づく各ペアのサンプル間のｂｒａｙ－Ｃｕｒｔｉｓの非類似度を算出した。以降の分析は、Ｒ（バージョン３．５．１）およびＪｕｐｙｔｅｒＬａｂ（バージョン０．３４．９）を用いて実施した。

【0164】

技術ノイズの推定
ＢａｒＢＩＱで測定した試料Ｍａ^ｄｉｓｔの技術反復中のｃＯＴＵのノイズは，Ｐｏｉｓｓｏｎ分布から得た模擬ノイズとｃＯＴＵ技術ノイズを比較することにより、主としてサンプリングノイズから確認した。技術反復における異なる検出総細胞数からのバイアスを排除するために、Ｒパッケージｖｅｇａｎにおける機能希薄化を用いたサブサンプリングにより、各反復の細胞数を、反復中の最小総細胞数に対して標準化した。ｃＯＴＵのノイズは、ＣＶ^２によって定量され、ここで、ＣＶは、３回の技術反復におけるｃＯＴＵの正規化細胞数に基づいて計算された係数の変動を示す^{５２，５３}。各ｃＯＴＵについてのシミュレートされたＰｏｉｓｓｏｎノイズを、試料中の与えられたｃＯＴＵの平均細胞数であるＰｏｉｓｓｏｎ分布からランダムに生成された３つの数字（３つの技術反復を模倣するため）に基づいて計算し、２つのシミュレーション（１および２）を行った。次に、各ｃＯＴＵについてＣＶ^２の理論平均補正後残差を算出した^{５２，５３}。
Ｒ_ｍｃ ＝ｌｏｇ_１０ （ＣＶ^２）－ｌｏｇ_１０ （ＣＶ_{Ｐｏｉｓｓｏｎ} ^２）；
ここで、ＣＶ_{Ｐｏｉｓｓｏｎ}は、ポアソン分布に基づく所定のｃＯＴＵに対する理論的ＣＶである。試料Ｍａ^ｄｉｓｔの全Ｒ^ｍｃの分布はシミュレーションの分布と一致しており、ＢａｒＢＩＱ測定の技術的ノイズは主としてサンプリングによるものであることが示唆された（図７ｃ，ｄ）。

【0165】

マウス依存性ＣＶの信頼区間推定
各ｃＯＴＵについて、シミュレーションにより、３匹のマウス（Ｍａ^ｄｉｓｔ１、Ｍｂ^ｄｉｓｔ、Ｍｃ^ｄｉｓｔまたはＭａ^ｐｒｏｘ１、Ｍｂ^ｐｒｏｘ、Ｍｃ^ｐｒｏｘ）の遠位または近位位置における細胞濃度のＣＶの９５％信頼区間を推定した。シミュレーションプロセスを１，０００回繰り返し、各時間について、所与のｃＯＴＵについて３つのシミュレートされた細胞濃度からＣＶを得た。各シミュレートされた濃度は、Ｐｏｉｓｓｏｎ分布から生成されたランダム数により得られ、その平均は、試料中の所与のｃＯＴＵのシークエンス決定細胞数（すなわち、Ｍａ^ｄｉｓｔ１、Ｍｂ^ｄｉｓｔ、Ｍｃ^ｄｉｓｔまたはＭａ^ｐｒｏｘ１、Ｍｂ^ｐｒｏｘ、Ｍｃ^ｐｒｏｘの１つ）であり、その後、この試料の推定総濃度を用いて正規化した。この推定総濃度は、その平均がこのサンプルの測定叢濃度であり、その標準偏差が平均の１０．１％である正規分布からランダムに生成した（１０．１％は、反復フィルタリングに関する５つの独立した実験の中で、平均によって標準化された最大標準偏差（１０．１％）であった（図１８））。各ＣＶの９５％信頼区間は、１，０００回のシミュレーションＣＶの分布から得た。

【0166】

修正細菌ネットワーク
階層的クラスタリングを、統計パッケージ中の機能ｈｃｌｕｓｔにより実施した（パッケージｐｈｅａｔｍａｐを用いてヒートマップを描いた）。クラスタリングに用いられる距離は、１－最小（│ｒ’│）［ｒ’∈（ｒ － ＯＣＩ， ｒ＋ＯＣＩ）］と定義され、ここで、ＯＣＩは各ｒの９０％片側信頼区間を意味する。階層的クラスタリングの系統樹を完全連結法により得た。具体的には、含まれるすべてのｃＯＴＵ間のピアソンの相関係数ｒを求めた。その後、ある微生物と他の微生物の距離を上記式に基づいて決定し、距離に基づいてｃＯＴＵをクラスタリングした。クラスタリング後の枝内の可能なｃＯＴＵペアの距離は枝の高さより低かった。各ｒのＯＣＩは、シミュレーションにより得た。シミュレーションプロセスを１，０００回繰り返し、各時間について、各ｃＯＴＵの細胞濃度を、試料Ｍａ^ｄｉｓｔ１、Ｍａ^ｐｒｏｘ１、Ｍｂ^ｄｉｓｔ、Ｍｂ^ｐｒｏｘ、Ｍｃ^ｄｉｓｔおよびＭｃ^ｐｒｏｘ（このプロセスは、上記のＣＶ信頼区間のシミュレーションと同じである）についてランダムに生成し、各ｃＯＴＵ対についてピアソンのｒを計算した。その後、ＯＣＩは、１，０００回のシミュレーションされたシミュレータの分布から得られた。

【0167】

閾値０．６により得られた各ＳＣＢＧについてのｃＯＴＵのネットワークを、パッケージのｉｇｒａｐｈを用いてｆｏｒｃｅ－ｄｉｒｅｃｔｅｄ ｌａｙｏｕｔ^３９により可視化し、ネットワークにおけるノード（すなわち、ｃＯＴＵ）のレイアウトを、０．９より大きいｒを用いて描き、ｃＯＴＵ間のすべてのｒを、パッケージＲＣｏｌｏｒＢｒｅｗｅｒを用いて色勾配により線で示した。

【0168】

ＳＣＢＧの各可能なペアの間のＳＣＢＧ間相関に基づくＳＣＢＧのネットワーク（Ｒ_{ｉｎｔｅｒ}）を、パッケージｉｇｒａｐｈを用いてｆｏｒｃｅ－ｄｉｒｅｃｔｅｄ ｌａｙｏｕｔにより可視化した。ＳＣＢＧのレイアウトは０．７より大きいＲ_{ｉｎｔｅｒ}に基づいて決定し、ＳＣＢＧ間のすべてのＲ_{ｉｎｔｅｒ}はパッケージＲＣｏｌｏｒＢｒｅｗｅｒを用いて色勾配により線で描いた。Ｒ_{ｉｎｔｅｒ}の平均の比較のためのＫｒｕｓｋａｌ－Ｗａｌｌｉｓ検定を、Ｒパッケージｓｔａｔｓにおける機能Ｋｒｕｓｋａｌ．ｔｅｓｔを用いて実施した。

【0169】

実施例１：細菌叢に含まれる細菌の単一細胞インデックス付与（ｉｎｄｅｘｉｎｇ）と単一ＲＮＡバーコード付与（Ｂａｒｃｏｒｄｉｎｇ）、ならびに配列解読による細胞単位と分子数のカウント

【0170】

微生物相と宿主の相互作用は宿主のホメオスタシスおよび多くの疾患^{１３－１６}と関連している。微生物相と宿主の相互作用のメカニズムをさらにかつ統合的に理解するために、微生物相を研究するだけでなく、微生物相と宿主の両方についてメタボロミクスおよび／またはトランスクリプトミクスのような他の分析と微生物相の組成分析を結びつけることが重要である^５。この目的のためには、一般に使用可能な単位、例えば、重量あたりの細胞数および／または体積あたりの分子数に基づく濃度測定が必要である。しかし、現行の技術^６－８では、細胞レベルでの微生物相組成の測定は困難であった。さらに、微生物叢は多数の細菌種からの膨大な数の細菌からなる^１７。そのため、高い分類学的分解能を有するハイスループット細胞定量法が望まれている。

【0171】

次世代シークエンシング技術を用いた１６Ｓ ｒＲＮＡ遺伝子増幅物シークエンシングに基づいたハイスループット法が、いく年にもわたる所与の細胞集団における細菌の多様性の研究に貢献した^{２２，２３}。しかし、従来の方法は、精製バルク細菌ゲノムから１６Ｓ ｒＲＮＡ遺伝子を増幅し、増幅した分子の数を測定するので、基本的には以下の限界がある。
１）異なる種はゲノム上に１６Ｓ ｒＲＮＡ遺伝子の異なるコピー数を有し、大部分の種のコピー数は不明であるため、細胞数を測定し、異なる種の細胞数を比較するのは困難である；
２）１６Ｓ ｒＲＮＡ配列の同定は、シークエンスおよび増幅エラーのために正確ではなく、低分類学的分解能をもたらす。
実際、シークエンスエラーは、分子バーコード^{２４－２６}を用いて修正されたが、主に配列増幅時に生じるキメラ生成に由来する増幅エラーは、未だ十分に除去できていない^２７。

【0172】

これらの従来の方法の限界を克服するために、正確な１６Ｓ ｒＲＮＡ遺伝子同定、ＢａｒＢＩＱ（図１ａ、図２９）を伴う細胞定量法を開発した。最初に試料を緩衝液中で調製し、ボルテックスにより塊を破壊した。次に、細菌試料を細胞バーコード^２５、^２６、プライマーおよびＤＮＡ増幅用試薬を含む溶液と混合し、１００μｍサイズの液滴に封入した。液滴の約４％が単一細胞と単一細胞バーコード（すなわち、ＤＮＡ分子）の両方を有するように、バーコードおよび細菌の濃度、ならびにポアソン分布に基づいてそれらの比を調整した。その後のシークエンスのために、増幅されたバーコードおよびシークエンスアダプターを、シングルステップ増幅によって液滴中の増幅された１６Ｓ ｒＲＮＡ遺伝子（Ｖ３－Ｖ４領域、約４５０塩基）と連結した^２８（図５）。増幅後、液滴を切断し、ライブラリー（連結したアンプリコン）を精製し、高スループットシーケンサーＭｉＳｅｑを用いて個々の増幅分子の細胞バーコードおよび１６Ｓ ｒＲＮＡ配列の両方をシークエンスした。われわれは、バーコード（すなわち、細胞）の各配列タイプについてシークエンスされた分子（すなわち、リード）を分析し、その１６Ｓ ｒＲＮＡ配列に基づいて各細胞のタイプを同定し、各細胞タイプについて細胞数をカウントした（図６）。この解析は、ゲノム上に複数の１６Ｓ ｒＲＮＡシークエンス型を有する細菌に対しても機能した。なぜなら、同じ細胞バーコードが、同一細胞由来の増幅された複数の１６Ｓ ｒＲＮＡシークエンスに連結したからである。最終的に、液滴デジタルＰＣＲ（「ＢａｒＢ１Ｑ法」の項目参照）で測定した同一試料の総濃度を用いて、シークエンスした細胞数を正規化することにより、試料中の各細胞型の細胞濃度を得た。

【0173】

ＢａｒＢＩＱと従来法との本質的な相違点は、微生物相の組成を定義するための単位である。従来法では、単位は操作上の分類単位（ｏｐｅｒａｔｉｏｎａｌ ｔａｘｏｎｏｍｉｃ ｕｎｉｔ， ＯＴＵ）であり、これはバルクサンプリングから得られた配列の同一性に基づいてクラスタリングすることによって基本的に得られる類似の１６Ｓ ｒＲＮＡ配列のグループを示す^３０。しかしながら、ＢａｒＢＩＱは、各バーコード化細胞から同定された１６Ｓ ｒＲＮＡ配列に基づいて分類される細胞型を使用する。本発明者らの細胞ベースの方法とＯＴＵを用いた従来の方法とを区別するために、ここで得られる分類単位を、「細胞ベースの操作上分類単位（ｃＯＴＵ）」と名付けた。

【0174】

最初に、ＢａｒＢＩＱが１０種類の培養ヒト腸内細菌株を含む模擬細胞集団に作用することを実証した（表１）。ＢａｒＢＩＱにより同定されたＢａｒ配列の２対を含む模擬細胞集団由来の１６の配列（Ｂａｒ配列）はそれぞれ１つの塩基差を有することを見出した（図１ｂ）。

【0175】

１６のＢａｒ配列はすべて、培養した１０株のサンガーシークエンスにより同定された１６Ｓ ｒＲＮＡ配列（Ｓａｎ配列）の１つと同一であった（図１ｂ）。いくつかのＳａｎ配列はＢａｒＢＩＱによっては見いだされず、サンガーシークエンスによって１個または２個の細胞からのみ検出されたことに注目する。次に、細胞バーコードに基づいて１６のＢａｒ－配列から１０のｃＯＴＵを同定し、それぞれが１０の株のうちの１つに対応した（図１ｂ、図２０～２８）。これに対して、通常の方法で同一の模擬細胞集団を測定したところ、１２種類の代表的なＯＴＵ配列のうち２種類のみがＳａｎ配列の１つと同一であることがわかった（図１ｂ）。したがって、ＢａｒＢＩＱは１６Ｓ ｒＲＮＡ配列同定のための一塩基精度および分解能を有した。ここで実施した従来の方法では実行不可能であると結論した。

【0176】

次に、ＢａｒＢＩＱによって、模擬細胞集団における各ｃＯＴＵの濃度（［Ｃ］_{ＢａｒＢＩＱ}） （体積当り）を測定した。ＢａｒＢＩＱで測定された濃度は、顕微鏡画像で測定されたそれらの細胞濃度（［Ｃ］_顕微鏡、図１ｃ）と一致することを確認した。２つの測定値間のＰｅａｒｓｏｎ積率相関係数ｒ（Ｐｅａｒｓｏｎのｒ）は０．９８であった。平均比（［Ｃ］_{ＢａｒＢＩＱ}／［Ｃ］_顕微鏡）は０．８８であり、固定勾配１（Ｒ^２＝０．９５）を用いて対数目盛に当てはめることにより決定した。この結果は、ＢａｒＢＩＱが模擬細胞集団における各細菌（ｃＯＴＵ）の細胞濃度を正確に測定したことを示唆した。

【0177】

次に、我々はＢａｒＢＩＱをマウス盲腸由来の微生物叢に適用した。盲腸は、微生物発酵容器^３１として機能し、微生物叢疾患関連研究のためのサンプリング場所としてしばしば選択されている^{３２，３３}。我々は、最近報告^３４されたように、細胞外細菌ＤＮＡが腸内微生物叢の定量に影響する可能性があるため、盲腸検体から細胞外ＤＮＡを除去した。

【0178】

Ｃ５７ＢＬ６／Ｊ雄マウス（Ｍａ、Ｍｂ、およびＭｃ）を共収容した３匹の結腸－盲腸および小腸－盲腸の両関節由来の２つの位置（遠位（ｄｉｓｔ）および近位（ｐｒｏｘ）で微生物叢を調査した（図２ａ）。Ｍａの両部位（遠位部位に対してＭａ^ｄｉｓｔ１、Ｍａ^ｄｉｓｔ２、Ｍａ^ｄｉｓｔ３および近位部位に対して；Ｍａ^ｐｒｏｘ１、Ｍａ^ｐｒｏｘ２およびＭａ^ｐｒｏｘ３のそれぞれ３回の反復実験を実施）から採取したサンプルと、その他（Ｍｂ^ｄｉｓｔ、Ｍｂ^ｐｒｏｘ、Ｍｃ^ｄｉｓｔおよびＭｃ^ｐｒｏｘ）について測定を実施した。合計で１．３×１０^５個の細菌細胞を計数し、７３０個のＢａｒ配列を含む６０４個のｃＯＴＵを同定した。驚くべきことに、同定された２３０のＢａｒ配列（７３０のうち３２％）は、広く使用されている３種類の公的データベース（ＧｒｅｅｎＧｅｎｅ^１０、Ｒｉｂｏｓｏｍａｌ Ｄａｔａｂａｓｅ Ｐｒｏｊｅｃｔ^１１、Ｓｉｌｖａ^１２）にいまだ登録されておらず、最も近い登録配列との同一性は８６．９～９９．９％であった（図２ｂ）。ＢａｒＢＩＱは単一塩基精度を有することを模擬細胞集団により明らかにしたため、ＢａｒＢＩＱは未知の１６Ｓ ｒＲＮＡ配列を同定することができると結論した。

【0179】

次に、試料中の同定された各ｃＯＴＵについて、上述のように細胞濃度を定量した。最初に、同一サンプルの技術的反復が再現性が高く（Ｐｅａｒｓｏｎのｒ≧０．９８２、図１５）、濾過工程を含む定量化のためのノイズは主にサンプリングによるものであることを確認した（図７のｂ－ｄ）。次いで、全サンプルで検出された２４０ｃＯＴＵを用いてサンプル間の細胞濃度を比較した。我々は、サンプルの各対（すなわち、異なる位置および／または異なるマウス）の１０～９７ｃＯＴＵ（図２ｄ、ｅ、および図１６ａ～１６ｃ）の濃度が異なることを見出した。その差はサンプリングノイズよりも大きく、倍率変化は２より大きかった（図２ｄ、ｅおよび図１６の信頼区間から外れ、倍率変化が２倍以上のドット；すなわち異なる濃度を有すると考えられるｃＯＴＵ）が、他の１４３～２３０（２４０ｃＯＴＵの６０～９６％）の濃度は一貫していた（図３ｄ、ｅおよび図１６ａ～ｃの信頼区間の範囲か、または倍率変化が２倍未満であるドット）。例えば、健康マウスと疾患マウスの間の差異を記述したいときは、一貫した細菌の同定が不可欠であろう。

【0180】

各サンプルペアの全体的な差を定量化するため、Ｂｒａｙ－Ｃｕｒｔｉｓの非類似性（存在量に基づくβの多様性）^３５分析を２４０ ｃＯＴＵの細胞濃度に基づいて実施した（図２ｆ）。上記の観察結果と一致して、異なるサンプル間の相違（図２ｆの異なる部位のサンプル間および異なるマウスからのサンプル間）は、Ｍａからの反復実験（図２ｆのＭａ^ｄｉｓｔおよびＭａ^ｐｒｏｘ参照）間の相違よりも有意に大きかった。さらに、近位位置または異なる位置（図２ｄの丸記号）での異なるマウスからのサンプル間の異質性は、同じマウスからの異なる位置、または異なるマウスからの遠位位置（図２ｄの三角記号）でのそれらより高かった。これらの結果は、微生物叢の細胞ベースの包括的な相違を定量的に分析し、同じマウスまたは異なるマウスの遠位と近位の位置からの微生物叢が全体的に異なることを示唆した。

【0181】

さらに、２４０ ｃＯＴＵの各々について、位置依存性（採取場所の位置（すなわち、遠位と近位）による依存性）およびマウス依存性の濃度差の両方を検討した。まず、Ｍａの各位置について３回の反復実験により、同じマウスの異なる位置間のｃＯＴＵ濃度を統計的に比較した。１３ ｃＯＴＵ（２４０ ｃＯＴＵのうち５％）は有意に異なり（ＦＤＲ＜０．０５かつ倍率変化＞２）（図８）、最大倍率変化（反復実験の平均濃度に基づく）は４．１であった。
次に、各部位のｃＯＴＵに対する変動係数（ＣＶ、３匹の細胞濃度の標準偏差をそれらの平均値で除したＣＶ）を算出することにより、３匹のマウスにおける細胞濃度の一貫性を定量した（図３ａ）。著者らは、各ｃＯＴＵについて遠位と近位の位置の間のＣＶを比較し（図９ａおよび９ｂ）、大部分のｃＯＴＵのＣＶは、シミュレートされた信頼区間に基づいて異ならないことを見出した。興味深いことに、同一属のｃＯＴＵｓの一貫性（すなわちＣＶ）は、しばしば異なっていた（図９ａおよび９ｂ）（公的データベースからの情報が限られているため、分類はＲＤＰ分類^３６を用いて各ｃＯＴＵのＢａｒ－ｓｅｑｕｅｎｃｅにより予測されたことに注目）。例えば、Ｃｌｏｓｔｒｉｄｉｕｍ ＸＩＶａ属のｃＯＴＵのＣＶは、両方の場所で０．０５～１．７０の範囲で変化した（図３ｂ）。興味深いことに、この属のある種は、盲腸^３１の主要な機能であるブチレート^３７のような短鎖脂肪酸を生成することが報告されている。この知見は、細菌の生理学的役割の更なる理解のためには、属レベルより細かいレベルの細胞の定量、特に、ｃＯＴＵレベルの細胞の定量が必要であることを示唆した^{３７，３８}。

【0182】

細胞間の関係を理解するために、各ｃＯＴＵ対の相関に基づいて細菌ネットワークを探索した。ヒト状態の遷移と関連する相関細菌ネットワークは何年かにわたって示されている。しかし、これまでのネットワーク分析は、基本的に、属レベルまたはより高いレベルでのＯＴＵ、すなわちｃＯＴＵではなくＯＴＵに基づいて行われた。本実施例では、６つの試料（Ｍａ^ｄｉｓｔ１、Ｍａ^ｐｒｏｘ１、Ｍｂ^ｄｉｓｔ、Ｍｂ^ｐｒｏｘ、Ｍｃ^ｄｉｓｔ、およびＭｃ^ｐｒｏｘ）からのそれらの測定濃度を用いて、一般的に検出される２９６ ｃＯＴＵの各対について対数目盛でピアソンのｒを計算することによって、細胞濃度に基づく相関を明らかにした（図４ａ）。大部分のｃＯＴＵと高い相関を示すｃＯＴＵは見いだされなかったが、いくつかのｃＯＴＵは他のいくつかのｃＯＴＵと強く相関した（図１０）。

【0183】

そこで、われわれは、可能性のあるすべてのｃＯＴＵ対の距離に基づいて階層的クラスタリングを行い、│ｒ│ｓを使用して、全てのｃＯＴＵが強く相関する細菌群（ｓｔｒｏｎｇｌｙ ｃｏｒｒｅｌａｔｅｄ ｂａｃｔｅｒｉａｌ ｇｒｏｕｐｓ； ＳＣＢＧｓ）を見出した（図４ｂ）。正および負の相関ｃＯＴＵの両方を含むＳＣＢＧを同定するために、ｒではなく│ｒ│を使用した、これは、ここで定義されたＳＣＢＧが「関係グループ」であることを意味する。ＳＣＢＧｓの同定の信頼性を確保するために、シミュレートされた誤差を│ｒ│ｓの計算に考慮した。図４ｂの系統樹上の点線として示される閾値０．６を使用したところ、計３１のＳＣＢＧが発見された。閾値より低く、３個以上のｃＯＴＵを含む枝をＳＣＢＧと定義した（図４ｂおよび図１７）。得られたＳＣＢＧの特徴を特徴付けた。ＳＣＢＧにおけるｃＯＴＵの数は３～１９個と様々であり（図１１ｃ）、半数以上（１６／３１）が正および負の相関したｃＯＴＵを含み（図４ｃおよび図１２ａ～ｆ）、各ＳＣＢＧにおけるｃＯＴＵの平均存在量が広く分布し（図４ｃおよび図１１ｄ）、ＳＣＢＧにおける存在量の最も高い差は２３０倍であり、ＳＣＢＧ１２（図１１ｄ）で認められ、（図３ｃと同様に）各ｃＯＴＧの分類学的予測（図４ｃ）はそれぞれ異種であった（図４ｃおよび図１３）。これらの知見は、細菌が異なる分類群に属するか、その存在量が異なる場合でも、強い関係を持つ可能性があることを示唆した。ＳＣＢＧの数および各ＳＣＢＧ中のｃＯＴＵの数は、閾値の関数として変化し（図１１ａおよび１１ｂ）、ここで使用した閾値０．６は、ＳＣＢＧの数に対する移行点であるように思われる（図１１ａ）。ＳＣＢＧを定義するために異なる閾値を選択して、特定の特徴が所定の試料において見出されるようにしてもよい。

【0184】

全ネットワークレベルでの細菌微生物叢の特性を評価するために、ＳＣＢＧの可能な全ての対を用いてＳＣＢＧ間の相関を調べた。２つのＳＣＢＧの間のＳＣＢＧ間の相関Ｒ_{ｉｎｔｅｒ}は、１つの対の２つのｃＯＴＵが異なるＳＣＢＧから由来する全ての可能なｃＯＴＵ対について計算された│ｒ│ｓの平均として定義した（図４ｄ）。最初に、全てのＲ_{ｉｎｔｅｒ}が、各ＳＣＢＧ内の全ての可能なｃＯＴＵペアについて計算された│ｒ│ｓの平均として定義されるＳＣＢＧ内の相関、Ｒ_{ｉｎｎｅｒ}よりも基本的に低いことを確認した（図１４ａ）。次いで、各ＳＣＢＧは少数のＳＣＢＧと相対的に高い相関を示すことを見出した。最後に、各ＳＣＢＧについて、ＳＣＢＧと他の全てのＳＣＢＧの間のＲ_{ｉｎｔｅｒ}の平均を計算し、全ての３１のＳＣＢＧの平均は、Ｋｒｕｓｋａｌ－Ｗａｌｌｉｓ検定を用いて有意差がないことを見出した（カイ二乗＝３０、ｄｆ＝２９、ｐ値＝０．４１）（図１４ｂ）。これらの知見は、ネットワーク全体におけるすべてのＳＣＢＧの平均的特徴は明確でないことを示唆した。全細菌ネットワークの分析は、疾患モデルの細菌ネットワークを健常マウスのそれと比較することにより、例えば、微生物叢における疾患関連ランドマークの特徴を見出すために重要であると思われる。

【0185】

更に、以下の検討を行った。
実験１． 大腸試料の細分化
大腸を細分化し、かつ細分化されたそれぞれの断片の位置情報を識別可能な状態で分析することを試みた。具体的には、マウスを屠殺した直後に大腸全体を摘出して直線状になるように広げ、それぞれの大腸固形内容物の位置関係を撮影記録した。一つの大腸固形内容物を腸壁に包まれた状態のまま滅菌済みハサミと滅菌済みピンセットを用いて取り出し、ブレインスライサ（室町機械、ＭＫ－ＲＣ－０１）の穴の中央に盲腸側を左側にして置いた（図３９のパネルａ）。このとき、盲腸側（Ａ）と肛門側（Ｅ）が識別できるように、ブレインスライサに予め印をつけた。

【0186】

次に、オートクレーブ滅菌後に５０℃で保温しておいた３％アガロース（ナカライテスク、０１１５７－９５）含有１×ＴＡＥ（ナカライテスク、３２６６６－８１）を静かに注ぎ入れ（図３９のパネルｂ）、－２０℃で３０分間静置することにより、大腸内容物をアガロースゲルに包埋した（図３９のパネルｃ）。ブレインスライサを－２０℃の冷凍庫から取り出し、大腸固形内容物の中央より１ｍｍほど左側に位置するブレインスライサの溝、およびその２つ右側に位置する溝に滅菌済みカミソリ刃（室町機械、ＴＣＢ－１００）を差し入れた（図３９のパネルｄ）。今回用いたブレインスライサの溝の幅は１ｍｍなので、この操作により大腸内容物から厚さ２ｍｍの中央部分（以下、Ｃエリア）が分割されたことになる。その後、大腸固形内容物の左端および右端から２ｍｍ以上３ｍｍ以下の内側に位置する溝にカミソリ刃をそれぞれ差し入れ（図３９のパネルｅ）、盲腸側先端部（以下、Ａエリア）および肛門側先端部（以下、Ｅエリア）に分けた。Ｂエリア（ＡエリアとＣエリアの間）およびＤエリア（ＣエリアとＥエリアの間）の幅は、大腸固形内容物全体の長さにより変動した。また、腸試料によっては内容物全体の長さが短く、ＢエリアまたはＤエリアのどちらか一方が欠落したものもあった。最後に、大腸固形内容物の左端および右端から１ｍｍ以上外側に位置する溝にカミソリ刃をそれぞれ差し入れた。カミソリ刃をブレインスライサから抜くことにより、各エリアの大腸内容物が含まれる切片をカミソリ刃に付着した状態で取り出し、滅菌済みピンセットを用いて各エリアの大腸内容物をＤＮＡ Ｌｏｂｉｎｄ Ｔｕｂｅ（Ｅｐｐｅｎｄｏｒｆ、００３０１０８０５１）に入れた。Ｃエリアは、ブレインスライサから取り出した後に滅菌済みペトリ皿に置いて（図３９のパネルｆ）、１５ゲージの滅菌済み金属ニードル（武蔵エンジニアリング、ＳＮＡ－１５Ｇ－Ｂ）を用いて中心部分（以下、ＣＣエリア）と周辺部分（以下、ＣＯエリア）に分け、ＤＮＡ Ｌｏｂｉｎｄ Ｔｕｂｅ に採取した。上記の一連の操作を、他の大腸固形内容物に対しても繰り返した。

【0187】

その結果、大腸固形内容物をＡ、Ｂ、ＣＣ、ＣＯ、Ｄ、Ｅの各エリアに分割することができた（図３９のパネルｇ）（ただし、上記理由によりＢまたはＤが欠落することがあった）。

【0188】

実験２． バーコード配列の分配に対する細菌の有無の影響の実験科学的検討
細菌細胞を含む、または含まない条件において、等モル混合した４つの細胞バーコードテンプレート（以下、等モル混合細胞バーコード）の濃度を、ｄｄＰＣＲを用いて測定した。具体的には、まず、ＱＸ２００^ＴＭ ｄｄＰＣＲ^ＴＭ ＥｖａＧｒｅｅｎ Ｓｕｐｅｒｍｉｘ（ＢｉｏＲａｄ、＃１８６４０３４）、１μＭプライマー（ＮｏＢｉｏｔｉｎ－Ｌｉｎｋ－ｂａｒｃｏｄｅ－ＦおよびＰ５－ｉｎｄｅｘ－Ｒ１Ｐ－ｂａｒｃｏｄｅ－Ｒ）、０．１μＭ ｄＮＴＰ（Ｎｅｗ Ｅｎｇｌａｎｄ ＢｉｏＬａｂｓ、Ｎ０４４７）、Ｐｌａｔｉｎｕｍ Ｔａｑ ＤＮＡ Ｐｏｌｙｍｅｒａｓｅ（Ｔｈｅｒｍｏ Ｆｉｓｈｅｒ Ｓｃｉｅｎｔｉｆｉｃ、１０９６６０３４）、およびサンプル（等モル混合細胞バーコードとマウス盲腸より採取した細菌細胞、あるいは等モル混合細胞バーコードのみ）を３０μｌの容量で混合し、ＤＧ８カートリッジ（ＢｉｏＲａｄ、＃１８６４００８）に分注した。次に、Ｄｒｏｐｌｅｔ Ｇｅｎｅｒａｔｉｏｎ Ｏｉｌ ｆｏｒ ＥｖａＧｒｅｅｎ（ＢｉｏＲａｄ、＃１８６４００６）およびＤｒｏｐｌｅｔ Ｇｅｎｅｒａｔｏｒ（ＢｉｏＲａｄ、＃１８６４００２ＪＡ）を用いて混合溶液を液滴に封入した。

【0189】

ｄｄＰＣＲ は以下のステップにより実施した。
第一ステージ；９５℃５分
第二ステージ；９５℃４５秒と６０℃１５０秒の繰り返しを６サイクル
第三ステージ；９５℃２５秒と６０℃８０秒の繰り返しを３９サイクル
第四ステージ；４℃５分と９０℃５分
その後、ＱＸ２００ Ｄｒｏｐｌｅｔ Ｒｅａｄｅｒ （ＢｉｏＲａｄ、＃１８６４００３ＪＡ）によりバーコードの濃度を測定した。

【0190】

結果、細菌細胞を含む、または含まない条件においてバーコード濃度の測定値に有意差はなく(図４０参照)、細菌細胞の有無はバーコードの液滴への分配比率に影響しないことが示された。

【0191】

実験３． ｄｄＰＣＲにおけるサイクル数の変更実験
ＢａｒＢＩＱ法のシークエンスライブラリー作製のためのＰＣＲサイクル数が、液滴に含まれる細菌細胞の１６Ｓ ｒＲＮＡ配列を増幅させるのに十分であることを確かめた。具体的には、まず、ＱＸ２００^ＴＭ ｄｄＰＣＲ^ＴＭ ＥｖａＧｒｅｅｎ Ｓｕｐｅｒｍｉｘ、１μＭプライマー（Ｆ１－ＦｗおよびＦ１－Ｒｖ）、０．１μＭ ｄＮＴＰ、およびサンプル（マウス盲腸より採取した細菌細胞）を３０μｌの容量で混合し、ＤＧ８カートリッジに分注した。
次に、Ｄｒｏｐｌｅｔ Ｇｅｎｅｒａｔｉｏｎ Ｏｉｌ ｆｏｒ ＥｖａＧｒｅｅｎおよび Ｄｒｏｐｌｅｔ Ｇｅｎｅｒａｔｏｒを用いて、混合溶液を液滴に封入した。ｄｄＰＣＲは、第三ステージを除き、上記実験２．と同様のサイクル条件で実施した。第三ステージは、サイクル数を０、１０、２０、３０、３９、あるいは４９に変えた。その後、液滴の蛍光強度をＱＸ２００ Ｄｒｏｐｌｅｔ Ｒｅａｄｅｒにより測定し、ソフトウェア ＱｕａｎｔａＳｏｆｔ（ＢｉｏＲａｄ、＃１８６４０１１ＪＡ）による強度の二峰性分布の谷である閾値に基づいて陽性および陰性液滴を決定した。

【0192】

結果、第三ステージのサイクル数が３０以上の条件下では、陽性液滴と陰性液滴の強度分布が明瞭に分離し（図４１のパネルａ）、液滴全体に占める陽性液滴の割合は１４％程度で一定となり（図４１のパネルｂ）、ＢａｒＢＩＱ法のシークエンスライブラリー作製のためのＰＣＲサイクル数が、液滴に含まれる細菌細胞の１６Ｓ ｒＲＮＡ配列を増幅させるのに十分であると示された。

【0193】

実験４． ｄｄＰＣＲのステップ時間の変更実験
ＢａｒＢＩＱ法のシークエンスライブラリー作製のための初期変性の時間が、液滴に含まれる細菌細胞の１６Ｓ ｒＲＮＡ配列を増幅させるのに十分であることを確かめた。具体的には、まず、ＱＸ２００^ＴＭ ｄｄＰＣＲ^ＴＭ ＥｖａＧｒｅｅｎ Ｓｕｐｅｒｍｉｘ、１μＭプライマー （Ｆ１－ＦｗおよびＦ１－Ｒｖ）、０．１μＭ ｄＮＴＰ、およびサンプル（マウス盲腸より採取した細菌細胞）を３０μｌの容量で混合し、ＤＧ８カートリッジに分注した。次に、Ｄｒｏｐｌｅｔ Ｇｅｎｅｒａｔｉｏｎ Ｏｉｌ ｆｏｒ ＥｖａＧｒｅｅｎおよびＤｒｏｐｌｅｔ Ｇｅｎｅｒａｔｏｒ を用いて、混合溶液を液滴に封入した。ｄｄＰＣＲは、第一ステージを除き、上記実験２．と同様のサイクル条件で実施した。第一ステージは、時間を０、５、あるいは１０分間に変えた。その後、液滴の蛍光強度をＱＸ２００ Ｄｒｏｐｌｅｔ Ｒｅａｄｅｒにより測定し、ソフトウェア ＱｕａｎｔａＳｏｆｔによる強度の二峰性分布の谷である閾値に基づいて陽性および陰性液滴を決定した。

【0194】

結果、第一ステージの時間を変えても、液滴全体に占める陽性液滴の割合は変化せず（図４２）、ＢａｒＢＩＱ法のシークエンスライブラリー作製のための初期変性の時間が、液滴に含まれる細菌細胞の１６Ｓ ｒＲＮＡ配列を増幅させるのに十分であると考えられた。

【0195】

現在、１６Ｓ ｒＲＮＡ遺伝子の増幅産物ベースの微生物叢の研究では、絶対定量^３、正確な測定^４０、完全遺伝子シークエンス^４１、および細菌－細菌相互作用を考慮する傾向にある^４２。しかし、これらは依然として細胞の定量とは関連していなかった。我々が知る限り、ＢａｒＢＩＱは、高スループット様式で細胞レベルでの細菌微生物相の高分類学的分解組成の定量を可能にする最初の方法である。さらに、ＢａｒＢＩＱによる未知の１６Ｓ ＲＮＡ配列の一塩基精度によるデータベースフリーの同定は、他の研究にとって有用であると思われる。例えば、新たに見出された細菌の局在を知りたいときは、ＢａｒＢＩＱにより同定された１６Ｓ ｒＲＮＡ配列を用いて蛍光プローブを設計することによって、ＦＩＳＨ（蛍光ｉｎ ｓｉｔｕハイブリダイゼーション）を行うことができる。

【0196】

最近、メタゲノミクス、トランスクリプトミクス、プロテオミクス、メタボロミクスのような異なるメタオミクスデータセットを統合し、これらのデータセットを用いたさらなる計算モデリングが、微生物相機能のメカニズム研究のための有望な方向として提案されている^５。このアプローチでは、細菌細胞は、明らかに異なるメタオミクスデータセットを統合するだけでなく、その機能のための基本単位であるので、微生物相は細胞レベルで定義されるべきである。ＢａｒＢＩＱによって提供される微生物叢の細胞ベースおよび分類群に依存しない解析は、微生物相の研究を現在の連合研究から必要な機構研究^４４に移行する。

【図1】

【図2】

【図3】

【図4】

【図5】

【図6】

【図7】

【図8】

【図9a】

【図9b】

【図10】

【図11】

【図12a】

【図12b】

【図12c】

【図12d】

【図12e】

【図12f】

【図13】

【図14】

【図15】

【図16a】

【図16b】

【図16c】

【図17】

【図18】

【図19】

【図20】

【図21】

【図22】

【図23】

【図24】

【図25】

【図26】

【図27】

【図28】

【図29】

【図30】

【図31】

【図32】

【図33】

【図34】

【図35】

【図36】

【図37】

【図38】

【図39】

【図40】

【図41】

【図42】

【配列表】

0007423101000001.app