異なるシーケンシング深さでのOtuのメンバーシップの変更(OTU不安定性)-微生物の多様性の分析のための無視されたが重要な特性
不安定なOtuによって作成された問題を説明するために、我々は同じデータセット(カナダ土壌データセット)と同じクラスタリング法(完全なリンケージクラスタリング、CLクラスタリングと呼ばれる)を用いて重複しない希薄化曲線を再現したRoeschらによって採用された。 （図1a）。 生のシーケンスを四つのシーケンスの深さでランダムにサブサンプ(20%, 40%, 60%, そして、入力配列の80％）それぞれの30回の複製を使用する。 次に、完全リンケージ（CL）クラスタリングを使用して、各サブサンプルをクラスター化し（すべてのクラスタリング方法の定義は追加のファイル1にあります）、各サンプリング深さに対して希薄化曲線を生成しました。 ＣＬクラスタリングの場合，より大きなサブサンプルによって生成される希薄化曲線は，より小さなサブサンプルによって生成される曲線よりも急である。
希薄化曲線を生成する際の目標の一つは、補間をサポートすることであり、完全なデータセットから希薄化曲線を作成する場合、その曲線を使用して、合計よりも少ない数のシーケンスで観測される種の数を決定したいと考えています。 たとえば、完全なデータセットから作成された希薄化曲線から補間すると、完全なデータセットから30,000個のシーケンスをランダムに選択すると、約4,500個の種があると推定されます（図1aのポイントA）。 しかし、重複しない希薄化曲線が補間に与える問題は、完全なデータセットの80%サブサンプルから30,000個のシーケンスをランダムにサブサンプリングした場合、これらの30,000個のシーケンスによって表されるのは4,200個だけであると推定されるということです(図1aのポイントB)。 このシナリオは、サンプルごとに少数の配列のみが収集された場合に本質的に当てはまり、希薄化曲線の予想される挙動と矛盾する現象である。
図1aに示すように、希薄化曲線の重複がないことは、実際にはOTUクラスタリング法の不安定性によって引き起こされることが観察されました。 つまり、シーケンスが割り当てられているクラスターは、クラスター化されているシーケンスの数によって影響を受ける可能性があります。 この仮説の図は図1bに示されています。s1とS2の2つの配列のみが類似しているしきい値（バーでリンクすることによって示されます）内で観察された場合、それらは単一のOTU（OTU1）にクラスター化されます。 次に、さらに3つのシーケンスS3、S4、およびS5を追加しますが、これはS1またはS2にリンクできますが、いくつかのペアワイズ距離がしきい値を超えます（これらのペアはバーでリンクされていません）。 CLの定義により、単一のOTUに割り当てられたすべての系列のペアワイズ距離は、距離閾値内に適合しなければならず、これにより、S1とS2がOTU2とOTU3に分離されることが可能になる可能性がある。 OTU1は、このシーケンスの深さで消え、その配列はOTU不安定性の問題を示す、二つの異なるOtuに再割り当てされます。 理論的には、より多くのシーケンスを追加すると、CLアルゴリズムを使用するときに既存のOtuが分割される傾向があります。 その結果、より大きなデータセットとより小さなデータセットでクラスター化されると、同じシーケンスがより多くのOtuにグループ化されます。 これにより、より大きなサンプルから得られる希薄化曲線が急峻になり、それがより高いアルファ多様性を有するという結論が得られる。 したがって、ＣＬから生じる希薄化曲線は、配列決定の深さにより敏感である。 この効果は弱いが、いくつかのケースでは、より小さなサンプルサイズに基づいているシーケンスの数の収集は、プラトーに到達する希薄化曲線を生成す 単にサンプリング深度の増加または減少のために異なるOtuに割り当てられている個人のこの現象は明らかに問題があります。 伝統的な（マクロスケールの）生態学に基づく類似の状況は、固定された領域内の異なる数の鳥を数えると、個々の鳥が種として一緒にグループ化する再定義 しかし、上記の不安定性は、伝統的な生態学の場合のように、新規種の時折の同定によるものではない。 対照的に、OTUメンバーシップへのこれらの変更は、Otu間で再割り当てされるシーケンスの大部分内で体系的に発生します。
不安定なOtuが生物学的解釈に及ぼす影響をさらに調査するために、次に序付けを使用してベータ多様性を検討しました。 主座標解析（PCoA）を使用して、我々は完全なデータセットの60％を含むサブサンプルを使用して完全なデータセットに対して微生物群集を比較した。 このサブサンプリングを30回繰り返して、反復を作成しました。 次に、CLクラスタリングを使用して、すべてのサブサンプルと完全なデータセットをクラスター化し、代表的なOTUシーケンス（各OTUで最も豊富なシーケンスと定義） その後、サンプルは、完全なデータセットのクラスタリングに起因する30個の複製の希少性を含む、サンプルあたり30,000個のシーケンスを含むように無作為に希少化されました。 希薄化後、すべてのサンプルには同じ数の配列が含まれていたため、それらの間の唯一の違いは最初にクラスター化された配列の数でした。 PCoAは、これらのサンプルが最初にクラスター化された配列の数に従って分離されていることを示し、OTUの不安定性が同じサンプルを異なる組成を有す 20％、40％、および80％のサブサンプルを完全なデータセットと比較したときにも同様の結果が観察されました（追加ファイル2：図S1）。 さらに、１ ２ ５Ｏｔｕ（偽発見率（ＦＤＲ）補正後）および２ ６Ｏｔｕ（Ｂｏｎｆｅｒｒｏｎｉ補正後）は、Ｍａｎｎ−Ｗｈｉｔｎｅｙ Ｕ試験を用いて、これら２つの群間で有意に異なることが決定された。 また、不安定なOtuが分類学的組成の計算に及ぼす影響をテストし、その効果は非常に限られていることがわかりました（追加ファイル3：図S2および追 これは、これらのOtuは、CLを使用してより多くの配列が追加されると変化しているにもかかわらず、系統発生的近接性の結果として同じ分類群に割り当てられているためである（不安定なOtuで系統発生指標を使用するPCoAの許容度を詳述するセクションでも後述する）。
代替の階層的および貪欲なクラスタリング方法は、不安定なOtuを生成します
OTUメンバーシップを決定するために使用されるすべての階層的方法は、Otu CLクラスタリングでは、距離しきい値内に収まるように、1つのOTU内のすべてのシーケンス間のペアワイズ距離が必要です。 単一リンケージ(SL)クラスタリングでは、距離しきい値内に収まるように、あるOTU内の任意のペアのシーケンス間のペアワイズ距離が必要です。 Average linkage(AL)クラスタリングでは、距離しきい値内に収まるように、1つのOTU内のすべてのシーケンス間の平均ペアワイズ距離が必要です。 SLクラスタリング(図2a)を使用する場合に予想されるように、Otuは、より多くのシーケンスが追加されるとマージされる傾向があります。これは、CLで観察される分割問題とは反対です。 したがって、SLを使用して作成された希薄化曲線は、サブサンプルサイズが大きくなるにつれて急勾配が少なくなります（図2b）。 ベータ多様性は、Otuの不安定なSLクラスタリングによっても影響されます(図2c)。 例えば、１ ６ ７Ｏｔｕ（ＦＤＲ補正後）および３ ６Ｏｔｕ（Ｂｏｎｆｅｒｒｏｎｉ補正後）は、６ ０％サブサンプルおよび完全データセットの両方にわたって差動的に表されることが決定された。
ブレイ-カーティス距離に基づく不安定な単一リンケージ（SL）クラスタリング、希薄化曲線、およびPCoAの基礎となる原則。 (a)二つのサンプリング深度での不安定なSLクラスタリングの基礎となる原則。 白い円は、小サブサンプルと大サブサンプルの両方に含まれていた個々のシーケンスを示し、暗い円は、大サブサンプルにのみ追加されたシーケンスを示し 線は、しきい値以下の距離を持つシーケンスのペアを示し、したがって、単一のOTUにリンクすることができます。 赤または青の大きな円は、それぞれ小サブサンプルと大サブサンプルのOTUsを示します。 (b,d)SL(b)と平均リンケージ(AL)(d)クラスタリングで生成された希薄化曲線。 （ｃ，ｅ）Ｂｒａｙ−Ｃｕｒｔｉｓ距離に基づくＰｃｏａ、ＳＬ（Ｃ）およびＡＬ（Ｅ）を用いて完全なデータセットと６ ０％のサブサンプルを比較する。 すべてのサブサンプルは、この分析に含まれるサンプルあたり30,000配列に希薄化されました。
OTU分割とOTUマージの両方が発生する可能性があるため、平均リンケージによって生成される不安定性はより複雑です。 これらの相反する効果は、OTUカウントのより微妙な違いにつながり、ALで作成された希薄化曲線は異なる深さで重複します（図2d）。 さらに、発生する多数のOTU分割およびマージイベントのために、AL Otu自体は不安定です（追加ファイル5：図S3）。 さらに、これらの不安定なOtuはベータ多様性に影響を与えるにもかかわらず（Adonis、R=0.16、P=0.001）、PCoAの主要な分離は、サンプルサイズ以外の要因によって引き起こ; たとえば、シーケンスの入力順序または異なるサブサンプル内の特定のキーシーケンスの有無に起因する違いを含めることができます（図2e）。 この観測は、入力シーケンスの順序に対するALの感度に起因する可能性があり、異なるクラスタリングパターンが生じる可能性があります。 ＡＬを使用する場合、８ ０ ４Ｏｔｕ（ＦＤＲ補正後）および５Ｏｔｕ（Ｂｏｎｆｅｒｒｏｎｉ補正後）は、２つのサンプリング深度にわたって差動的に表された。
Usearchで実装されているような貪欲なクラスタリングは、CL、SL、ALよりも計算効率が高いde novoクラスタリング方法としてよく使用されています。 欲張りクラスタリングを使用する場合、シーケンスは、そのOTU内でクラスター化される単一のOTU重心の距離しきい値内にある必要があります。 さらに、シーケンスは定義された順序で処理され、各クエリシーケンスは既存のOTUに割り当てられるか、新しいOTUの重心として割り当てられます。 一つのクエリシーケンスが複数の既存のOTU重心の距離しきい値内にある場合、最も近い重心（ここでは距離ベースの貪欲クラスタリング（DGC）と呼ばれる）また しかし、私たちは、最も一般的に採用されているものに焦点を限定することを選択しました。 本研究では，貪欲クラスタリングの方法としてＵＳＥＡＲＣＨを評価した（ｕｐａｒｓｅのクラスタリングアルゴリズムはＵＳＥＡＲＣＨで使用されているものと同じであるため，ＵＰＡＲＳＥは評価しなかった）。
OTU不安定性は貪欲なクラスタリング法の問題でもあり、いくつかのソースから発生します。 第一に、重心の選択は、配列が処理される順序に大きく依存する。 したがって、サンプルのサイズが変更されると、配列の順序も変更され得る。 第二に、DGCを使用する場合、サンプルのサイズが増加したときに重心の選択が安定したままであっても、追加された配列は新しい重心になり、既存のOtuからメンバーを引き付けることができる（これは一般的にAGCでは起こらない）。 たとえば、s10、S11、およびS12が、S10を重心としてOTU7を形成するとします（図3a、b）。 後続の配列決定において、別の配列Ｓ１ ３が追加された場合、より大きなサンプルの処理順序は、Ｓ１ ０、Ｓ１ ３、Ｓ１ １、およびＳ１ ２となり得る。 この場合、S10は依然として重心になりますが、S13も重心になります。 そして、ｓ１ ３は、２つの間の距離がＳ１ １とＳ１ ０との間の距離よりも小さいので、Ｓ１ １を募集する。 DGCでは、S11はS10ではなくS13でクラスタリングされ、元のOTU7はOTU8とOTU9に分割されます（図3a）。 AGCでは、S11はまだS10とクラスタ化され、元のOTU7は元の構造を保持します（図3b）。
不安定距離ベースの貪欲クラスタリング（DGC）と豊富ベースの貪欲クラスタリング（AGC）、希薄化曲線、およびブレイ-カーティス距離に基づくPCoAの基礎となる原則。 (a,b)二つのサンプリング深さで不安定なDGC(a)とAGC(b)の基礎となる原理。 白い円は、小サブサンプルと大サブサンプルの両方に含まれていた個々のシーケンスを示し、暗い円は、大サブサンプルにのみ追加されたシーケンスを示し 黄色の点はOTU重心を示しています。 線は、しきい値以下の距離を持つシーケンスのペアを示し、したがって、単一のOTUにリンクすることができます。 赤または青の大きな円は、それぞれ小サブサンプルと大サブサンプルのOTUsを示します。 （ｃ，ｄ）Ｄｇｃ（ｃ）とＡＧＣ（ｄ）で五つの異なる深さで生成された希薄化曲線。 （ｅ、ｆ）Ｂｒａｙ−Ｃｕｒｔｉｓ距離に基づくＰｃｏａ、Ａ ＧＣ（Ｅ）およびＤＧＣ（ｆ）を用いて完全なデータセットと６ ０％のサブサンプルを比較する。 すべてのサブサンプルは、この分析に含まれるサンプルあたり30,000配列に希薄化されました。
我々は、不安定なOtuによって生成される効果を分析するために、アルファ希薄化曲線とベータ多様性PCoA上の貪欲クラスタリングを使用しました。 前述したように、DGCとAGCはどちらも重心の変更性に苦しんでいます（この効果はOTUの分割またはマージに偏っていません）。 その結果、DGCとCLクラスタリングは同様の曲線を生成し、サブサンプルのサイズが大きくなるにつれて急峻になりました（図3c）。 対照的に、AGCは深さの影響を受けない重複した曲線を生成しました（図3d）。 しかし，ＡＬクラスタリングと同様に，これはＯｔｕが安定していることを意味するものではなく，異なるサブサンプリング深度で同様の数の（おそらく異なる）Ｏｔｕが得られたことだけを意味する。 ベータ多様性のDGCおよびAGC効果推定で生成された不安定なOtu（図3e、f）。 ＡＧＣの場合、３ ９ ２Ｏｔｕ（ＦＤＲ補正後）および１ ４Ｏｔｕ（Ｂｏｎｆｅｒｒｏｎｉ補正後）は、２つの深さにわたって差動的に表されると決定され、ＤＧＣの場合、これらの数は、それぞれ３ ７ ０および１ ５であ
これらの不安定な方法の違いを定量化するために、不安定な配列と不安定なOtuの割合を比較しました（図4a、b、追加ファイル6：表S1）。 CLは不安定な配列（約22％）の最も高い割合を生産し、AL（13％）とAGC（12％）はSL（15％）とDGC（14％）よりもわずかに優れていた。 これらの結果は、代替データセットの使用を比較するときに常に一貫しているわけではありませんでした(追加ファイル7:図S6);しかし、AGCは一般に、他のde novo法と比較して最高のパフォーマンスを示しました。 不安定なOtuについては、CLとDGCが不安定なOtuの割合が最も高く、周波数が10以上の重心を持つOtuの約60%がそれぞれの方法で不安定であることが観察された(>90%は、追加のファイル7:図S6に示すように、特定のデータセットを分析すると不安定であることが判明した)。 ALおよびSLはCLまたはDGCのいずれかよりも安定であるが、それでも少なくとも30回観察されている重心のための30％以上のOTU不安定性をもたらした。 ＡＧＣは，特に高度に豊富な重心を持つＯｔｕに対して，最も安定なｄｅｎｏｖｏ法であることが分かった。
不安定な配列の割合、不安定なOtuの割合、および各メソッドのMCC値。 (a)方法によって作成された不安定な配列の割合。 不安定なシーケンスは、60%サブサンプルでは1つの重心にクラスター化されていますが、100%(完全)データセットでは別の重心にクラスター化されているシーケンスと定義されます。 (b)方法およびクラスタ重心の頻度によって作成された不安定なOtuの割合(閉じた参照と逆複製の値はゼロであり、この図には含まれていません)。 OTUが60%と100%のデータセットで同一であった場合(60%サブサンプルに存在しないシーケンスは含まない)、それは安定として定義されます。 （ｃ）各方法のＭＣＣ値。 より高い値は、より高い安定性に対応する。
安定したOtuを生成するde novoクラスタリング方法の1つは、逆複製または同一で長さが等しいシーケンスのクラスタリングです（追加ファイル8：図S4A）。 Closed-reference OTUクラスタリングと同様に、クラスタリングはクラスター化されるシーケンスコレクションの構成によって影響を受けないため、すべてのOtuは異な その結果、デレプリケーションを使用して生成された希薄化曲線は、異なる深さにわたって重複しています（追加ファイル8: 図S4B)、ベータダイバーシチはサブサンプルのサイズの影響を受けません(追加ファイル8:図S4C)。 さらに、単一のＯＴＵが２つの群間で有意に異なると決定されるわけではない。 デレプリケーションは、シーケンスエラーに起因する偽のOtuを識別することに対して非常に脆弱であることに注意することが重要です。 ビニングOtuの安定性のために、それはまた、不安定なOtuが（シーケンスエラーではなく）重複していない希薄化曲線の主な原因であることを示す、異なる深さ さらに，デノボクラスタリング法が類似性しきい値が増加するにつれてデノボクラスタリング法がデレプリケーションに類似するようになるので，デレプリケーション法の安定性は，クラスタリングのための高い類似性しきい値が不安定なＯｔｕの発生を減少させることを示唆している。 実際には、デレプリケーションクラスタリングは、下流を使用するために計算上高価であるOtuの高い数をもたらします。 したがって、現代のデータセットサイズは、逆複製されただけのシーケンスを扱うことを妨げています。 今後の方法では，ＯＴＵ不安定性の問題を管理するために逆解法に基づくアプローチを使用する可能性がある。 もう1つの極端な例は、すべてのシーケンスを1つのOTUにクラスタリングすることですが、そのOTUは絶対に安定しています。 それにもかかわらず、デレプリケーションとは異なり、OTUsは、α多様性、β多様性、および分類学的組成などのさらなる分析に利用することができる。 さらに、すべてのシーケンスを一つのOTUにクラスタリングすることは”クラスタリング”とはほとんど呼ばれず、下流の分析には完全に役に立たない。
参照ベースの方法は、不安定なOtuの問題を最小限に抑える
すべての不安定なクラスタリング方法が共通している一つの特徴は、クラスター定義が入力シーケンスに依存していることです。 クローズ参照OTUクラスタリングでは、事前に計算された重心の参照データセット(Greengenesデータベースなど)に対して読み取りがクラスタリングされ、クラスタリング中に新しい重心が作成されないため、完全に安定したOtuが生成されます(図5a)。 その結果、閉じた参照クラスタリングに基づくアルファとベータの多様性の推定は、サンプルのサイズの影響を受けず（図5b、c）、Otuは2つの深さの間で有意に異なると判断されることはありません。 安定したOtuを生成することに加えて、閉じた参照クラスタリングは、他のいくつかの便利な機能を提供します。 第一に、参照配列の名前は、任意に割り当てられた名前を使用するのではなく、普遍的なＯＴＵ識別子として使用することができ、したがって、研究間のＯｔｕの直接比較を容易にする。 第二に、参照データセットが完全長のマーカー遺伝子で構成されている場合、異なるマーカー遺伝子領域からの配列読み取りを一緒にクラスター化するこ 最後に、閉じた参照クラスタリングは、大規模なデータセットのOTUクラスタリングを並列化できます。 Closed-reference OTUクラスタリングの主な制限は、任意の参照重心との類似性しきい値の外側にある読み取りが破棄され、データベース内で既に表されているOtuのみが”カナダの土壌データセットの処理では、配列の約14％が参照配列に一致することができなかったため、クラスタリング後に破棄されました。 このようなクローズドリファレンスOTUクラスタリングの制限は、参照データセットに対する改善が予想されるため、特定の研究プロジェクト（例えば、腸ミクロバイオーム）に必要な対応する参照がより高度に開発されるようになるため、些細なものになる可能性があります。
ブレイ-カーティス距離に基づく安定した閉基準クラスタリング、希薄化曲線、およびPCoAの基礎となる原則。 (a)二つのサンプリング深さで安定した閉じた参照クラスタリングの基礎となる原則。 白い円は、小サブサンプルと大サブサンプルの両方に含まれていた個々のシーケンスを示し、暗い円は、大サブサンプルにのみ追加されたシーケンスを示し ダイヤモンドは参照配列を示します。 線は、しきい値以下の距離を持つシーケンスのペアを示し、したがって、単一のOTUにリンクすることができます。 赤または青の大きな円は、それぞれ小サブサンプルと大サブサンプルのOTUsを示します。 (b)五つの異なる深さで閉じた参照クラスタリングで生成された希薄化曲線。 (c)PCoAは、Bray-Curtis距離に基づいて、閉じた参照クラスタリングを使用して完全なデータセットと60％のサブサンプルを比較します。 すべてのサブサンプルは、この分析に含まれるサンプルあたり30,000配列に希薄化されました。
クローズ参照OTUクラスタリングの制限を克服するために、オープン参照OTUクラスタリングを使用することができます。 オープン参照クラスタリングは、クローズ参照クラスタリングと同じ方法で開始されますが、参照コレクションと一致しないシーケンスをde novoの方法でクラ 既存のde novoクラスタリング法は不安定なOtuを生成しますが、多くのシーケンスが最初に閉じた参照アプローチによってクラスター化されるため、オープン参照クラスタリングはそのような方法よりもはるかに安定している可能性があります。 我々は、de novoクラスタリングステップ（図4a、b、c）のAGCを使用してオープン参照クラスタリングにおけるOTUの安定性を評価し、それがde novo法のみを使用するよりもはるかに効果的な方法であることがわかった。 不安定なOtuの大部分は、参照一致のない低存在量配列（一般的にエラーが発生しやすいと考えられている配列のカテゴリ）でした。 オープンリファレンスOTUクラスタリングは重複希薄化曲線を生成し（追加ファイル9：図S5A）、オープンリファレンスOTUクラスタリングの不安定性がPCoA解析に影響を与えるにもかかわらず（追加ファイル9: ０ ３）は、２つのグループにわたって差動的に表されるＯｔｕの数（ＦＤＲ補正後の１ ０ ４ＯｔｕおよびＢｏｎｆｅｒｒｏｎｉ補正後の２Ｏｔｕ）と同様に、任意の他のｄｅ ｎｏｖｏ法単独の場合よりも低い。 我々は、不安定なシーケンスと不安定なOtuの割合に焦点を当て、追加のデータセット上の他のde novoメソッドとオープンリファレンスクラスタリング法を比較し、これらの結果は、環境の種類とシーケンシング技術にわたって一般的に一貫していることがわかった(追加ファイル7:図S6)。
Otuの不安定性を定量化することに加えて、MCCインデックスを使用して、完全なデータセットと60%サブセットのクラスタリングに基づいてシーケンスペアのクラスタリングがどのように変化したかを調査しました(図4b、追加ファイル6:表S2)。 二つの参照ベースの方法と逆複製クラスタリングは、このメトリックによって最も高い安定性を有し、AGCはde novoクラスタリング法の中で最も安定であることは明らかである(Kruskal-Wallis test,P<0.05)。 ＡＬはｍｃｃ値が最も低く，この方法を用いたときに多くのシーケンスペアのクラスタリングが変化したことを示した。 あるいは，ＳＬはＡ ｌおよびＣＬを含むｄｅｎｏｖｏ法のほとんどよりも高いＭＣＣ値を生成した。 それにもかかわらず、SLのMCC値が高い理由の一部は、そのFP値が0に等しいことです（小さいサブサンプルで分離されたシーケンスは、大きなサブサンプ したがって、OTUマージに関する深刻な問題のために、SLははるかに安定した方法とはみなされるべきではありません。
系統発生ベータダイバーシティメトリクスはOTU不安定性の影響を最小限に抑える
すべてのOtuが互いに等しくないと考えられる非系統発生メトリクスとは異なり、UniFracのような系統発生メトリクスは、サンプル間の距離を計算する際にOtu間の系統発生関係を考慮に入れている。 不安定なOTUクラスタリング法は、通常は進化的に密接に関連しているOtu間のシーケンスを移動するため、サンプル間の計算された距離は、一般的に非系統 UniFrac距離に基づいてCL、SL、AL、AGC、およびDGCを使用して、ベータ多様性に対する不安定なOtuの効果を再分析しました（追加ファイル10：図S7）。 結果は、CL、AGC、およびDGCの不安定なOtuは、最小限のUniFrac距離を使用してベータ多様性に影響を与えることを示し、配列がこれらの不安定な方法と密接に関連したOtuの間で変化しているとき、系統発生メトリックは、その不安定性に対してより耐性があるという仮説を確認している。 それにもかかわらず、SLクラスタリングでは、遠くに関連するOtuを最終的に単一のOTUに結合することができるため、UniFrac距離を使用している場合でもベータダイバーシティが影響を受ける可能性があります。 ALでは、主要な分離は、非系統発生指標と同様に、異なるクラスタリングパターンによって引き起こされます。