ミオシン重鎖を超えたヒト骨格筋線維の不均一性

nature.com をご利用いただきありがとうございます。お使いのブラウザのバージョンは、CSS のサポートが制限されています。最適な表示のために、最新バージョンのブラウザをご利用いただくか、Internet Explorer の互換表示モードを無効にすることをお勧めします。また、継続的なサポートを確保するため、このサイトではスタイルシートと JavaScript は使用しません。
骨格筋は、主に筋原線維から構成される不均一な組織であり、ヒトでは筋原線維は通常、1つの「遅筋」（タイプ1）と2つの「速筋」（タイプ2Aおよび2X）の3種類に分類されます。しかし、従来の筋原線維タイプ間およびタイプ内の不均一性については、まだ十分に理解されていません。私たちは、ヒト大腿外側広筋からそれぞれ1050個と1038個の筋原線維にトランスクリプトミクスとプロテオミクスのアプローチを適用しました。プロテオミクス研究には男性が含まれ、トランスクリプトミクス研究には男性10人と女性2人が含まれました。ミオシン重鎖アイソフォームに加えて、代謝タンパク質、リボソームタンパク質、および細胞間結合タンパク質が、筋原線維間の多次元的な変動の原因であることを明らかにしました。さらに、遅筋線維と速筋線維のクラスターが特定されたにもかかわらず、私たちのデータは、タイプ2X線維が他の速筋線維と表現型的に区別できないことを示唆しています。さらに、ミオシン重鎖に基づく分類だけでは、ネマリンミオパチーにおける筋線維の表現型を十分に説明することはできません。全体として、我々のデータは、ミオシン重鎖アイソフォーム以外にも変動要因が存在する、多次元的な筋線維の不均一性を示唆しています。
細胞の不均一性は、すべての生物システムに固有の特徴であり、細胞が組織や細胞のさまざまなニーズを満たすように特殊化することを可能にする。1 骨格筋線維の不均一性に関する従来の見解は、運動ニューロンが運動単位内の線維タイプを規定し、線維タイプ（すなわち、ヒトではタイプ1、タイプ2A、およびタイプ2X）はミオシン重鎖（MYH）アイソフォームの特性によって決定されるというものであった。2 これは当初、それらのpH ATPase不安定性3,4に基づいており、後にMYHの分子発現5に基づいていた。しかし、さまざまな割合で複数のMYHを共発現する「混合」線維の同定とその後の受容により、骨格筋線維は、明確な線維タイプではなく、連続体としてますます見なされるようになっている。6 それにもかかわらず、この分野では依然として筋線維分類の主要な分類因子としてMYHに大きく依存しており、この見解は、MYH発現プロファイルと線維タイプの範囲がヒトとは異なる初期のげっ歯類研究の限界と重大なバイアスに影響されている可能性が高い。ヒトの場合も同様である。2 状況は、ヒトの骨格筋の種類によって線維の種類が異なるという事実によってさらに複雑になる。7 大腿外側広筋は、中間的な（したがって代表的な）MYH発現プロファイルを持つ混合筋である。7 さらに、サンプリングが容易であるため、ヒトで最もよく研​​究されている筋肉となっている。
したがって、強力な「オミクス」ツールを用いた骨格筋線維の多様性の偏りのない調査は重要であると同時に、骨格筋線維が多核であることなどから困難でもある。しかし、近年の様々な技術進歩により、トランスクリプトミクス8,9およびプロテオミクス10技術は感度において革命的な進歩を遂げ、単一線維解像度での骨格筋の解析が可能になった。その結果、単一線維の多様性や、萎縮刺激および加齢に対するそれらの応答の特性評価において大きな進歩が遂げられた11,12,13,14,15,16,17,18。重要なことに、これらの技術進歩は臨床応用が可能であり、疾患に関連する調節異常をより詳細かつ正確に特徴づけることができる。例えば、最も一般的な遺伝性筋疾患の1つであるネマリンミオパチー（MIM 605355およびMIM 161800）の病態生理は複雑で不明瞭です。19,20 したがって、骨格筋線維の調節異常をより詳細に特徴づけることで、この疾患の理解が大きく進展する可能性があります。
私たちは、ヒト生検標本から手作業で分離した単一の骨格筋線維のトランスクリプトーム解析およびプロテオーム解析法を開発し、それを数千本の線維に適用することで、ヒト骨格筋線維の細胞的不均一性を調査しました。この研究の過程で、筋線維のトランスクリプトーム解析およびプロテオーム解析による表現型解析の有効性を実証し、代謝タンパク質、リボソームタンパク質、および細胞間結合タンパク質が線維間変動の重要な原因であることを特定しました。さらに、このプロテオーム解析ワークフローを用いて、単一の骨格筋線維における線虫ミオパチーの臨床的意義を特徴づけ、MYHに基づく線維タイプに関係なく、非酸化型線維への協調的なシフトを明らかにしました。
ヒト骨格筋線維の不均一性を調査するために、単一の骨格筋線維のトランスクリプトームおよびプロテオーム解析を可能にする 2 つのワークフローを開発しました (図 1A および補足図 1A)。サンプルの保存と RNA およびタンパク質の完全性の維持から、各アプローチのスループットの最適化まで、いくつかの方法論的ステップを開発および最適化しました。トランスクリプトーム解析では、逆転写の最初のステップでサンプル固有の分子バーコードを挿入することにより、96 本の線維をプールして効率的な下流処理を行うことができました。従来の単一細胞アプローチと比較して、より深いシーケンス (線維あたり ± 100 万リード) により、トランスクリプトームデータがさらに充実しました。21 プロテオミクスでは、短いクロマトグラフィー勾配 (21 分) と timsTOF 質量分析計での DIA-PASEF データ取得を組み合わせて、高いスループットを維持しながらプロテオーム深度を最適化しました。 22,23 健康な骨格筋線維の不均一性を調査するために、14 人の健康な成人ドナーから採取した 1,050 本の線維のトランスクリプトームと、5 人の健康な成人ドナーから採取した 1,038 本の線維のプロテオームを特徴付けました (補足表 1)。この論文では、これらのデータセットをそれぞれ 1,000 本の線維のトランスクリプトームとプロテオームと呼びます。私たちのアプローチでは、1,000 本の線維のトランスクリプトームとプロテオームの解析で合計 27,237 個の転写産物と 2,983 個のタンパク質を検出しました (図 1A、補足データセット 1–2)。トランスクリプトームとプロテオームのデータセットを、検出された遺伝子が 1,000 個以上、かつ線維あたり 50% の有効値を持つようにフィルタリングした後、トランスクリプトームとプロテオームのそれぞれ 925 本の線維と 974 本の線維に対して、後続のバイオインフォマティクス解析を実行しました。フィルタリング後、平均して4257 ± 1557個の遺伝子と2015 ± 234個のタンパク質（平均 ± SD）が線維ごとに検出され、個体間の変動は限定的でした（補足図1B～C、補足データセット3～4）。ただし、被験者内の変動は参加者間でより顕著であり、これはおそらく長さや断面積の異なる線維間のRNA/タンパク質収量の違いによるものです。ほとんどのタンパク質（>2000）の変動係数は20%未満でした（補足図1D）。両方の方法により、筋収縮に重要な高発現シグネチャを持つ転写産物とタンパク質の広いダイナミックレンジを捉えることができました（例：ACTA1、MYH2、MYH7、TNNT1、TNNT3）（補足図1E～F）。特定された特徴のほとんどは、トランスクリプトームデータセットとプロテオームデータセットの間で共通しており（補足図1G）、これらの特徴の平均UMI/LFQ強度は、かなり良好な相関関係を示しました（r = 0.52）（補足図1H）。
トランスクリプトミクスおよびプロテオミクスのワークフロー（BioRender.comで作成）。BD MYH7、MYH2、MYH1のダイナミックレンジ曲線と、ファイバータイプの割り当てのために計算された閾値。E、F トランスクリプトミクスおよびプロテオミクスデータセットにおけるファイバー全体にわたるMYH発現の分布。G、H MYHベースのファイバータイプで色分けされたトランスクリプトミクスおよびプロテオミクスのUniform Diversity Approximation and Projection（UMAP）プロット。I、J トランスクリプトミクスおよびプロテオミクスデータセットにおけるMYH7、MYH2、MYH1の発現を示す特徴プロット。
我々は当初、オミクスデータセットにおけるMYH発現の高い感度とダイナミックレンジを活用する最適化されたアプローチを用いて、各線維にMYHベースの線維タイプを割り当てることに着手した。これまでの研究では、異なるMYHの発現の固定パーセンテージに基づいて、線維を純粋なタイプ1、タイプ2A、タイプ2X、または混合として分類するために任意の閾値が使用されてきた11,14,24。我々は、線維のタイプ分けに使用したMYH（MYH7、MYH2、MYH1、それぞれタイプ1、タイプ2A、タイプ2X線維に対応）によって各線維の発現をランク付けするという異なるアプローチを使用した。次に、得られた各曲線の下側の変曲点を数学的に計算し、それを閾値として使用して、各MYHについて線維を陽性（閾値以上）または陰性（閾値以下）として割り当てた（図1B～D）。これらのデータは、MYH7 (図 1B) と MYH2 (図 1C) が、タンパク質レベルと比較して RNA レベルでより明確なオン/オフ発現プロファイルを持っていることを示しています。実際、タンパク質レベルでは、MYH7 を発現していない線維はごくわずかで、MYH2 を 100% 発現している線維はありませんでした。次に、あらかじめ設定した発現閾値を使用して、各データセットのすべての線維に MYH ベースの線維タイプを割り当てました。たとえば、MYH7+/MYH2-/MYH1- 線維はタイプ 1 に割り当てられ、MYH7-/MYH2+/MYH1+ 線維は混合タイプ 2A/2X に割り当てられました (詳細については補足表 2 を参照)。すべての線維をプールすると、予想どおり、MYH ベースの線維タイプの分布は RNA (図 1E) とタンパク質 (図 1F) の両方で非常に類似していることが観察されましたが、MYH ベースの線維タイプの相対的な構成は個人によって異なりました (補足図 2A)。ほとんどの線維は純粋なタイプ1（34～35％）またはタイプ2A（36～38％）に分類されたが、かなりの数の混合型タイプ2A/2X線維も検出された（16～19％）。顕著な違いは、純粋なタイプ2X線維はRNAレベルでは検出できたが、タンパク質レベルでは検出できなかったことであり、これはMYHの急速な発現が少なくとも部分的には転写後レベルで制御されていることを示唆している。
我々は、抗体を用いたドットブロッティング法を用いて、プロテオミクスに基づくMYH線維型判定法の妥当性を検証したところ、両手法とも純粋なタイプ1線維とタイプ2A線維の同定において100%の一致率を達成した（補足図2B参照）。しかしながら、プロテオミクスに基づく手法は、より感度が高く、混合線維の同定や各線維における各MYH遺伝子の割合の定量化においてより効率的であった。これらのデータは、骨格筋線維型を特徴づけるために、客観的で高感度なオミクスに基づく手法を用いることの有効性を示している。
次に、トランスクリプトミクスとプロテオミクスによって提供される情報を組み合わせて、筋線維の完全なトランスクリプトームまたはプロテオームに基づいて客観的に分類しました。均一多様体近似および投影（UMAP）法を使用して次元を6つの主成分に削減し（補足図3A～B）、トランスクリプトーム（図1G）とプロテオーム（図1H）における筋線維の変動を可視化することができました。注目すべきは、トランスクリプトミクスデータセットとプロテオミクスデータセットのどちらにおいても、筋線維は参加者（補足図3C～D）またはテスト日（補足図3E）によってグループ化されておらず、骨格筋線維の被験者内変動が被験者間変動よりも大きいことを示唆していることです。UMAPプロットでは、「速い」筋線維と「遅い」筋線維を表す2つの明確なクラスターが現れました（図1G～H）。 MYH7+（遅筋）筋線維はUMAP1の正極に集積し、MYH2+およびMYH1+（速筋）筋線維はUMAP1の負極に集積した（図1I～J）。しかし、MYH発現に基づいて速筋線維タイプ（すなわち、タイプ2A、タイプ2X、または混合2A/2X）を区別することはできず、MYH1（図1I～J）またはACTN3やMYLK2（補足図4A～B）などの他の古典的な2X筋線維マーカーの発現は、トランスクリプトームまたはプロテオーム全体を考慮すると、異なる筋線維タイプを区別しないことが示唆された。さらに、MYH2およびMYH7と比較すると、MYH1と正の相関を示す転写産物またはタンパク質は少なく（補足図4C～H）、MYH1の発現量が筋線維のトランスクリプトーム/プロテオームを完全に反映しているわけではないことが示唆された。3つのMYHアイソフォームの混合発現をUMAPレベルで評価した場合も同様の結論が得られた（補足図4I～J）。したがって、MYHの定量化のみに基づいて転写産物レベルで2X線維を識別することは可能であるが、トランスクリプトームまたはプロテオーム全体を考慮すると、MYH1+線維は他の速筋線維と区別できない。
MYH以外の遅筋線維の不均一性の初期調査として、確立された4つの遅筋線維型特異的タンパク質、TPM3、TNNT1、MYL3、およびATP2A22を評価しました。遅筋線維サブタイプは、トランスクリプトミクス（補足図5A）とプロテオミクス（補足図5B）の両方で、MYH7と高い（ただし完全ではない）ピアソン相関を示しました。遅筋線維の約25%と33%は、それぞれトランスクリプトミクス（補足図5C）とプロテオミクス（補足図5D）のすべての遺伝子/タンパク質サブタイプによって純粋な遅筋線維として分類されませんでした。したがって、複数の遺伝子/タンパク質サブタイプに基づく遅筋線維の分類は、線維型特異的であることが知られているタンパク質であっても、さらなる複雑さをもたらします。これは、単一の遺伝子/タンパク質ファミリーのアイソフォームに基づく線維分類が、骨格筋線維の真の不均一性を適切に反映していない可能性があることを示唆しています。
ヒト骨格筋線維の表現型の多様性をオミクスモデル全体のスケールでさらに詳しく調べるために、主成分分析（PCA）を使用してデータのバイアスのない次元削減を実行しました（図2A）。UMAPプロットと同様に、参加者もテスト日もPCAレベルでの線維クラスタリングに影響を与えませんでした（補足図6A～C）。両方のデータセットで、MYHベースの線維タイプはPC2によって説明され、遅筋型1線維のクラスターと、速筋型2A、2X、および混合2A/2X線維を含む2番目のクラスターを示しました（図2A）。両方のデータセットで、これら2つのクラスターは少数の混合型1/2A線維によって接続されていました。予想どおり、主要なPCドライバーの過剰表現分析により、PC2は収縮および代謝シグネチャによって駆動されていることが確認されました（図2Bおよび補足図6D～E、補足データセット5～6）。全体として、MYHに基づく線維タイプは、PC2に沿った連続的な変動を説明するのに十分であることがわかった。ただし、高速クラスター内のトランスクリプトーム全体に分布しているいわゆる2X線維は例外である。
A. MYH に基づく線維タイプに応じて色分けされたトランスクリプトームおよびプロテオームデータセットの主成分分析 (PCA) プロット。 B. PC2 および PC1 における転写産物およびタンパク質ドライバーの濃縮分析。統計分析は clusterProfiler パッケージと Benjamini-Hochberg 調整済み p 値を使用して実行されました。 C、D. トランスクリプトームの細胞間接着遺伝子オントロジー (GO) 用語とプロテオームのコスタメア GO 用語に応じて色分けされた PCA プロット。矢印は転写産物およびタンパク質ドライバーとその方向を表します。 E、F. 遅筋線維/速筋線維タイプとは無関係な発現勾配を示す、臨床的に関連する特徴の均一多様体近似および投影 (UMAP) 特徴プロット。 G、H. トランスクリプトームおよびプロテオームにおける PC2 および PC1 ドライバー間の相関。
予想外にも、MYH に基づく筋線維タイプは、2 番目に高い変動度 (PC2) しか説明できず、MYH に基づく筋線維タイプ (PC1) とは無関係の他の生物学的因子が骨格筋線維の不均一性を制御する上で重要な役割を果たしていることが示唆された。PC1 の上位ドライバーの過剰表現分析により、PC1 の変動は主にトランスクリプトームの細胞間接着とリボソーム含有量、およびプロテオームのコスタメアとリボソームタンパク質によって決定されることが明らかになった (図 2B および補足図 6D–E、補足データセット 7)。骨格筋では、コスタメアは Z ディスクを筋形質膜に接続し、力の伝達とシグナル伝達に関与している。25 細胞間接着 (トランスクリプトーム、図 2C) およびコスタメア (プロテオーム、図 2D) の特徴を使用した注釈付き PCA プロットは、PC1 の強い左シフトを示し、これらの特徴が特定の線維に豊富に存在することを示している。
UMAPレベルでの筋線維クラスタリングのより詳細な解析により、ほとんどの特徴が筋線維サブクラスター特異的ではなく、筋線維タイプに依存しないMYHベースの発現勾配を示すことが明らかになった。この連続性は、CHCHD10（神経筋疾患）、SLIT3（筋萎縮）、CTDNEP1（筋疾患）など、病理学的状態に関連するいくつかの遺伝子で観察された（図2E）。この連続性は、神経疾患に関連するタンパク質（UGDH）、インスリンシグナル伝達（PHIP）、転写（HIST1H2AB）など、プロテオーム全体でも観察された（図2F）。これらのデータは総合的に、異なる筋線維間で、線維タイプに依存しない遅筋/速筋の不均一性の連続性を示している。
興味深いことに、PC2 のドライバー遺伝子はトランスクリプトームとプロテオームの相関が良好 (r = 0.663) (図 2G) であり、遅筋線維と速筋線維のタイプ、特に骨格筋線維の収縮特性と代謝特性が転写的に制御されていることを示唆している。しかし、PC1 のドライバー遺伝子はトランスクリプトームとプロテオームの相関が全くなく (r = -0.027) (図 2H)、遅筋線維と速筋線維のタイプとは無関係な変動は主に転写後制御されていることを示唆している。PC1 の変動は主にリボソーム遺伝子オントロジー用語で説明され、リボソームはタンパク質翻訳に積極的に関与し影響を与えることで細胞内で重要かつ特殊な役割を果たしていることから、31 次にこの予期せぬリボソームの不均一性を調査することにした。
まず、GOCC用語「細胞質リボソーム」に含まれるタンパク質の相対存在量に基づいて、プロテオミクス主成分分析プロットに色を付けました（図3A）。この用語はPC1の正の側に濃縮されているため、勾配は小さいですが、リボソームタンパク質はPC1の両方向の分割を駆動しています（図3A）。PC1の負の側に濃縮されたリボソームタンパク質には、RPL18、RPS18、およびRPS13が含まれていました（図3B）。一方、RPL31、RPL35、およびRPL38（図3C）は、PC1の正の側の主な駆動因子でした。興味深いことに、RPL38とRPS13は、他の組織と比較して骨格筋で高発現していました（補足図7A）。PC1におけるこれらの特徴的なリボソームシグネチャはトランスクリプトームでは観察されなかったため（補足図7B）、転写後調節が示唆されます。
A. プロテオーム全体にわたる細胞質リボソーム遺伝子オントロジー（GO）用語に従って色分けされた主成分分析（PCA）プロット。矢印は、PCAプロットにおけるタンパク質を介した変動の方向を示します。線の長さは、特定のタンパク質の主成分スコアに対応します。B、C. RPS13とRPL38のPCA特徴プロット。D. 細胞質リボソームタンパク質の非教師あり階層的クラスタリング分析。E. 骨格筋線維で異なる存在量を持つリボソームタンパク質を強調した80Sリボソームの構造モデル（PDB：4V6X）。F. mRNA出口チャネル付近に局在する、異なる化学量論比を持つリボソームタンパク質。
リボソームの不均一性と特殊化の概念は以前から提唱されており、異なるリボソームサブポピュレーション（リボソームの不均一性）の存在が、特定のmRNA転写プールの選択的翻訳34（リボソームの特殊化）を介して、異なる組織32および細胞33におけるタンパク質翻訳に直接影響を与える可能性がある。骨格筋線維で共発現するリボソームタンパク質のサブポピュレーションを特定するために、プロテオーム中のリボソームタンパク質の非教師あり階層的クラスタリング分析を行った（図3D、補足データセット8）。予想どおり、リボソームタンパク質はMYHに基づく線維タイプによってクラスター化されなかった。しかし、リボソームタンパク質の3つの異なるクラスターを特定した。最初のクラスター（ribosomal_cluster_1）はRPL38と共調節されており、正のPC1プロファイルを持つ線維で発現が増加している。2番目のクラスター（ribosomal_cluster_2）はRPS13と共調節されており、負のPC1プロファイルを持つ線維で増加している。 3 番目のクラスター (ribosomal_cluster_3) は骨格筋線維で協調的な差次的発現を示さず、「コア」骨格筋リボソームタンパク質とみなすことができます。リボソームクラスター 1 と 2 の両方には、以前に代替翻訳を調節することが示され (例: RPL10A、RPL38、RPS19、および RPS25)、機能的に発生に影響を与える (例: RPL10A、RPL38) リボソームタンパク質が含まれています。34,35,36,37,38 PCA の結果と一致して、これらのリボソームタンパク質の線維全体にわたる異質な表現も連続性を示しました (補足図 7C)。
リボソーム内の異種リボソームタンパク質の位置を可視化するために、ヒト80Sリボソームの構造モデル（Protein Data Bank: 4V6X）を使用しました（図3E）。異なるリボソームクラスターに属するリボソームタンパク質を単離した後、それらの位置は密接に一致しておらず、このアプローチではリボソームの特定の領域/フラクションを濃縮できなかったことが示唆されました。しかし興味深いことに、クラスター2の大きなサブユニットタンパク質の割合は、クラスター1および3よりも低くなりました（補足図7D）。骨格筋線維で化学量論が変化したタンパク質は主にリボソーム表面に局在していることが観察されました（図3E）。これは、異なるmRNA集団の内部リボソーム進入部位（IRES）要素と相互作用し、選択的翻訳を調整する能力と一致しています。 40, 41 さらに、骨格筋線維で化学量論が変化した多くのタンパク質は、mRNA出口トンネル（図3F）などの機能領域の近くに位置しており、特定のペプチドの翻訳伸長と停止を選択的に制御している。42 要約すると、我々のデータは、骨格筋リボソームタンパク質の化学量論が不均一であり、骨格筋線維間で差異が生じることを示唆している。
次に、速筋線維と遅筋線維の特徴を特定し、それらの転写調節のメカニズムを解明することに着手した。2 つのデータセット (図 1G–H および 4A–B) で UMAP によって定義された速筋線維と遅筋線維のクラスターを比較すると、トランスクリプトーム解析とプロテオーム解析により、それぞれ 1366 個と 804 個の発現量が異なる特徴が特定された (図 4A–B、補足データセット 9–12)。サルコメア (トロポミオシンやトロポニンなど)、興奮収縮連関 (SERCA アイソフォーム)、エネルギー代謝 (ALDOA や CKB など) に関連する特徴に予想通りの違いが見られた。さらに、タンパク質ユビキチン化を調節する転写産物やタンパク質 (USP54、SH3RF2、USP28、USP48 など) が速筋線維と遅筋線維で異なって発現していた (図 4A–B)。さらに、以前に子羊の筋線維タイプ間で発現が異なることが示され43、心筋でSERCA活性を高める44微生物タンパク質遺伝子RP11-451G4.2（DWORF）が、遅筋線維で有意に上方制御されていた（図4A）。同様に、個々の線維レベルでは、代謝関連の乳酸脱水素酵素アイソフォーム（LDHAおよびLDHB、図4Cおよび補足図8A）45,46などの既知のシグネチャ、およびこれまで知られていなかった線維タイプ特異的シグネチャ（IRX3、USP54、USP28、およびDPYSL3など）に有意な差が観察された（図4C）。トランスクリプトームデータセットとプロテオームデータセットの間には、発現が異なる特徴の有意な重複があり（補足図8B）、また、主にサルコメアの特徴のより顕著な発現の差によって引き起こされる倍率変化の相関があった（補足図8C）。特に注目すべきは、一部のシグネチャー（USP28、USP48、GOLGA4、AKAP13など）は、プロテオームレベルでのみ強い転写後調節を示し、遅筋線維/速筋線維タイプ特異的な発現プロファイルを示した（補足図8C）。
A および B ボルケーノ プロットは、図 1G ～ H の均一多様体近似および投影 (UMAP) プロットによって識別された遅いクラスターと速いクラスターを比較したものです。色の付いたドットは、FDR < 0.05 で有意に異なる転写産物またはタンパク質を表し、濃い色のドットは、log 変化 > 1 で有意に異なる転写産物またはタンパク質を表します。二元統計解析は、DESeq2 Wald 検定と Benjamini–Hochberg 調整済み p 値 (トランスクリプトミクス)、または Limma 線形モデル法と経験的ベイズ分析、それに続く多重比較のための Benjamini–Hochberg 調整 (プロテオミクス) を使用して実行されました。C 遅い繊維と速い繊維の間で選択された差次的発現遺伝子またはタンパク質のシグネチャ プロット。D 有意に差次的発現した転写産物とタンパク質の濃縮分析。重複する値は両方のデータセットで濃縮され、トランスクリプトームの値はトランスクリプトームのみで濃縮され、プロテオームの値はプロテオームのみで濃縮されます。統計解析は、Benjamini-Hochberg法で調整したp値を用いたclusterProfilerパッケージを使用して実施した。E. SCENIC由来の調節因子特異性スコアと線維タイプ間のmRNA発現差に基づいてSCENICによって同定された線維タイプ特異的転写因子。F. 遅筋線維と速筋線維間で発現差のある選択された転写因子のプロファイリング。
次に、差次的発現遺伝子とタンパク質の過剰発現解析を実施しました（図4D、補足データセット13）。2つのデータセット間で異なる特徴のパスウェイエンリッチメントにより、脂肪酸β酸化とケトン代謝プロセス（遅筋線維）、筋フィラメント/筋収縮（それぞれ速筋線維と遅筋線維）、炭水化物異化プロセス（速筋線維）など、予想された違いが明らかになりました。セリン/スレオニンタンパク質ホスファターゼ活性も速筋線維で上昇しており、これはグリコーゲン代謝を調節することが知られている調節および触媒ホスファターゼサブユニット（PPP3CB、PPP1R3D、およびPPP1R3A）などの特徴によって引き起こされています（47）（補足図8D～E）。速筋線維で豊富に見られた他の経路には、プロテオームにおけるプロセシング（P-）ボディ（YTHDF3、TRIM21、LSM2）（補足図8F）、転写後調節に関与する可能性のあるもの（48）、およびトランスクリプトームにおける転写因子活性（SREBF1、RXRG、RORA）（補足図8G）が含まれていた。遅筋線維では、酸化還元酵素活性（BDH1、DCXR、TXN2）（補足図8H）、アミド結合（CPTP、PFDN2、CRYAB）（補足図8I）、細胞外マトリックス（CTSD、ADAMTSL4、LAMC1）（補足図8J）、および受容体リガンド活性（FNDC5、SPX、NENF）（補足図8K）が豊富であった。
遅筋線維と速筋線維のタイプ特性の根底にある転写調節についてさらに理解を深めるため、SCENIC49 を使用した転写因子濃縮解析を実施した（補足データセット 14）。多くの転写因子が速筋線維と遅筋線維の間で有意に濃縮されていた（図 4E）。これには、以前に速筋線維の発達に関連付けられていた MAFA などの転写因子 50 や、これまで筋線維タイプ特異的遺伝子プログラムに関連付けられていなかったいくつかの転写因子が含まれていた。これらのうち、PITX1、EGR1、および MYF6 は速筋線維で最も濃縮された転写因子であった（図 4E）。対照的に、ZSCAN30 および EPAS1（HIF2A としても知られる）は遅筋線維で最も濃縮された転写因子であった（図 4E）。これと一致して、MAFA は速筋線維に対応する UMAP 領域でより高いレベルで発現していたが、EPAS1 は逆の発現パターンを示した（図 4F）。
既知のタンパク質コード遺伝子に加えて、ヒトの発達や疾患の調節に関与する可能性のある多数の非コードRNAバイオタイプが存在する。51, 52 トランスクリプトームデータセットでは、いくつかの非コードRNAが線維型特異性を示し（図5Aおよび補足データセット15）、その中には遅筋線維に非常に特異的で、ミトコンドリアミオパチー患者の筋肉で減少していることが報告されているLINC01405が含まれる。53 対照的に、lnc-ERCC5-5遺伝子（https://lncipedia.org/db/transcript/lnc-ERCC5-5:2）54に対応するRP11-255P5.3は速筋線維型特異性を示す。 LINC01405（https://tinyurl.com/x5k9wj3h）とRP11-255P5.3（https://tinyurl.com/29jmzder）はどちらも骨格筋特異性を示し（補足図9A～B）、1 Mbのゲノム近傍に既知の収縮遺伝子が存在しないことから、これらの遺伝子は近傍の収縮遺伝子を制御するのではなく、線維タイプの制御に特化した役割を果たしていることが示唆される。LINC01405とRP11-255P5.3の遅筋線維/速筋線維タイプ特異的発現プロファイルは、RNAscopeを用いて確認された（図5B～C）。
A. 非コードRNA転写産物は、遅筋線維と速筋線維で有意に調節されている。B. LINC01405とRP11-255P5.3の遅筋線維と速筋線維のタイプ特異性をそれぞれ示す代表的なRNAscope画像。スケールバー = 50 μm。C. RNAscopeによって決定された筋線維タイプ特異的非コードRNA発現の定量化（n = 3人の独立した個人からの生検、各個人内の速筋線維と遅筋線維を比較）。統計分析は、両側スチューデントt検定を使用して実行した。ボックスプロットは中央値と第1四分位数と第3四分位数を示し、ひげは最小値と最大値を指している。D. 新規微生物タンパク質同定ワークフロー（BioRender.comで作成）。 E. 微生物タンパク質 LINC01405_ORF408:17441:17358 は、遅筋線維に特異的に発現しています (n=5 個の独立した参加者からの生検サンプル、各参加者の速筋線維と遅筋線維を比較)。統計解析は、Limm 線形モデル法と経験的ベイズ法を組み合わせて実施し、続いて多重比較のための Benjamini-Hochberg 法と p 値の調整を行いました。ボックスプロットは中央値、第 1 四分位数、第 3 四分位数を示し、ひげは最大値/最小値を示しています。
最近の研究では、多くの推定非コード転写産物が転写された微生物タンパク質をコードしており、その一部が筋肉機能を調節していることが示されています。44, 55 潜在的な線維型特異性を持つ微生物タンパク質を特定するために、1000 線維トランスクリプトームデータセットで見つかった非コード転写産物 (n = 305) の配列を含むカスタム FASTA ファイルを使用して、1000 線維プロテオームデータセットを検索しました (図 5D)。22 種類の転写産物から 197 種類の微生物タンパク質を特定し、そのうち 71 種類は遅筋線維と速筋線維の間で異なる調節を受けていました (補足図 9C および補足データセット 16)。LINC01405 については、3 つの微生物タンパク質産物が特定され、そのうちの 1 つは転写産物と同様の遅筋線維特異性を示しました (図 5E および補足図 9D)。したがって、LINC01405 は遅筋線維に特異的な微生物タンパク質をコードする遺伝子であると特定しました。
我々は、個々の筋線維の大規模なプロテオミクス特性評価のための包括的なワークフローを開発し、健康な状態における線維の不均一性の調節因子を特定した。我々は、このワークフローをネマリンミオパチーが骨格筋線維の不均一性にどのように影響するかを理解するために適用した。ネマリンミオパチーは、筋力低下を引き起こす遺伝性筋疾患であり、罹患した子供では、呼吸困難、脊柱側弯症、四肢可動域の制限など、さまざまな合併症を呈する。19,20 一般的に、ネマリンミオパチーでは、アクチンα1（ACTA1）などの遺伝子の病原性変異により、遅筋線維の筋線維構成が優勢になるが、この効果は不均一である。注目すべき例外の1つは、速筋線維が優勢なトロポニンT1ネマリンミオパチー（TNNT1）である。したがって、ネマリンミオパチーで観察される骨格筋線維の機能障害の根底にある異質性をより深く理解することで、これらの疾患と筋線維の種類との複雑な関係を解明するのに役立つ可能性がある。
健常対照群（各群n=3）と比較して、ACTA1およびTNNT1遺伝子に変異を有するネマリンミオパチー患者から分離した筋線維は、顕著な筋線維萎縮またはジストロフィーを示した（図6A、補足表3）。利用可能な材料の量が限られているため、プロテオミクス解析には重大な技術的課題が生じた。それにもかかわらず、272本の骨格筋線維で2485個のタンパク質を検出することができた。線維あたり少なくとも1000個の定量化されたタンパク質をフィルタリングした後、250本の線維をバイオインフォマティクス解析に供した。フィルタリング後、線維あたり平均1573±359個のタンパク質が定量化された（補足図10A、補足データセット17～18）。注目すべきは、線維サイズが大幅に減少したにもかかわらず、ネマリンミオパチー患者サンプルのプロテオーム深度はわずかに減少しただけであったことである。さらに、独自のFASTAファイル（非コード転写産物を含む）を使用してこれらのデータを処理することで、ネマリンミオパチー患者の骨格筋線維から5つの微生物タンパク質を同定することができました（補足データセット19）。プロテオームのダイナミックレンジは有意に広く、コントロール群の総タンパク質は、以前の1000線維プロテオーム解析の結果とよく相関していました（補足図10B～C）。
A. ACTA1およびTNNT1ネマリンミオパチー（NM）におけるMYHに基づく線維萎縮またはジストロフィーおよび異なる線維タイプの優位性を示す顕微鏡画像。スケールバー = 100 μm。ACTA1およびTNNT1患者における染色再現性を確保するため、代表的な画像を選択する前に、3人の患者の生検組織を2～3回染色した（1症例あたり4切片）。B. MYHに基づく参加者の線維タイプの割合。C. ネマリンミオパチー患者および対照群の骨格筋線維の主成分分析（PCA）プロット。 D. 図 2 で解析した 1000 本の線維から決定した PCA プロット上に投影された、ネマリンミオパチー患者と対照群の骨格筋線維。たとえば、ACTA1 および TNNT1 ネマリンミオパチー患者と対照群の間、および ACTA1 および TNNT1 ネマリンミオパチー患者間の差異を比較したボルケーノプロット。色付きの円は π < 0.05 で有意差のあるタンパク質を示し、黒い点は FDR < 0.05 で有意差のあるタンパク質を示します。統計解析は Limma 線形モデル法と経験的ベイズ法を使用して実行され、その後、Benjamini-Hochberg 法を使用して多重比較の p 値調整が行われました。 H. 全プロテオームおよびタイプ 1 および 2A 線維における有意に発現差のあるタンパク質の濃縮解析。統計解析は clusterProfiler パッケージと Benjamini-Hochberg 調整済み p 値を使用して実行されました。 I、J. 細胞外マトリックスとミトコンドリア遺伝子オントロジー（GO）用語で色分けされた主成分分析（PCA）プロット。
ネマリンミオパチーは骨格筋におけるMYH発現筋線維の割合に影響を与える可能性があるため、19,20 まずネマリンミオパチー患者と対照群におけるMYH発現筋線維の種類を調べた。1000筋線維アッセイについて以前に説明したバイアスのない方法を使用して筋線維の種類を決定したが（補足図10D～E）、純粋な2X筋線維を再び特定できなかった（図6B）。ネマリンミオパチーが筋線維の種類に及ぼす影響は不均一であることが観察され、ACTA1変異を有する2人の患者ではタイプ1筋線維の割合が増加したが、TNNT1ネマリンミオパチーを有する2人の患者ではタイプ1筋線維の割合が減少した（図6B）。実際、ACTA1ネマリンミオパチーではMYH2および速筋トロポニンアイソフォーム（TNNC2、TNNI2、TNNT3）の発現が低下し、TNNT1ネマリンミオパチーではMYH7の発現が低下していた（補足図11A）。これは、ネマリンミオパチーにおける異質な筋線維タイプのスイッチングに関する以前の報告と一致している。19,20 我々は免疫組織化学によってこれらの結果を確認し、ACTA1ネマリンミオパチーの患者ではタイプ1筋線維が優勢であったのに対し、TNNT1ネマリンミオパチーの患者では逆のパターンが見られた（図6A）。
単一線維プロテオームレベルでは、ACTA1およびTNNT1ネマリンミオパチー患者の骨格筋線維は、大部分のコントロール線維とクラスターを形成し、TNNT1ネマリンミオパチー線維が一般的に最も重度に影響を受けていました（図6C）。これは、各患者の擬似膨張線維の主成分分析（PCA）プロットをプロットすると特に顕著で、TNNT1ネマリンミオパチー患者2と3はコントロールサンプルから最も離れているように見えました（補足図11B、補足データセット20）。ミオパチー患者の線維が健康な線維とどのように比較されるかをよりよく理解するために、健康な成人参加者の1,000本の線維のプロテオーム分析から得られた詳細な情報を使用しました。ミオパチーデータセット（ACTA1およびTNNT1ネマリンミオパチー患者とコントロール）の線維を、1,000本の線維のプロテオーム分析から得られたPCAプロットに投影しました（図6D）。対照線維におけるPC2に沿ったMYH線維タイプの分布は、1000線維プロテオーム解析から得られた線維分布と類似していた。しかし、ネマリンミオパチー患者の線維のほとんどは、本来のMYH線維タイプに関わらず、PC2を下方にシフトし、健康な速筋線維と重なっていた。したがって、ACTA1ネマリンミオパチー患者はMYHベースの方法を用いて定量化した場合、タイプ1線維へのシフトを示したが、ACTA1ネマリンミオパチーとTNNT1ネマリンミオパチーの両方において、骨格筋線維プロテオームは速筋線維の方へシフトしていた。
次に、各患者群を健常対照群と直接比較し、ACTA1およびTNNT1ネマリンミオパチーでそれぞれ256および552の差次的発現タンパク質を同定した（図6E～Gおよび補足図11C、補足データセット21）。遺伝子エンリッチメント解析により、ミトコンドリアタンパク質の協調的な減少が明らかになった（図6H～I、補足データセット22）。驚くべきことに、ACTA1およびTNNT1ネマリンミオパチーにおける線維タイプの優位性の違いにもかかわらず、この減少はMYHベースの線維タイプとは完全に独立していた（図6Hおよび補足図11D～I、補足データセット23）。ACTA1またはTNNT1ネマリンミオパチーでは、3つの微生物タンパク質も調節されていた。これらのマイクロプロテインの 2 つ、ENSG00000215483_TR14_ORF67 (別名 LINC00598 または Lnc-FOXO1) と ENSG00000229425_TR25_ORF40 (lnc-NRIP1-2) は、タイプ 1 筋線維でのみ異なる存在量を示しました。ENSG00000215483_TR14_ORF67 は、以前に細胞周期制御に関与することが報告されています。56 一方、ENSG00000232046_TR1_ORF437 (LINC01798 に対応) は、ACTA1 ネマリンミオパチーのタイプ 1 およびタイプ 2A 筋線維の両方で、健常対照と比較して増加していました (補足図 12A、補足データセット 24)。対照的に、リボソームタンパク質はネマリンミオパチーの影響をほとんど受けなかったが、ACTA1ネマリンミオパチーではRPS17の発現が低下していた（図6E）。
エンリッチメント解析では、ACTA1およびTNNT1ネマリンミオパチーで免疫系プロセスのアップレギュレーションが明らかになった一方、TNNT1ネマリンミオパチーでは細胞接着も増加していた（図6H）。これらの細胞外因子のエンリッチメントは、PC1およびPC2のPCAを負の方向（すなわち、最も影響を受けた線維の方向）にシフトさせる細胞外マトリックスタンパク質によって反映された（図6J）。両患者群とも、アネキシン（ANXA1、ANXA2、ANXA5）57,58およびそれらの相互作用タンパク質S100A1159など、免疫応答および筋形質膜修復機構に関与する細胞外タンパク質の発現増加を示した（補足図12B～C）。このプロセスは、以前に筋ジストロフィーで増強されることが報告されているが60、我々の知る限り、ネマリンミオパチーと関連付けられたことはない。この分子機構の正常な機能は、損傷後の筋形質膜の修復と、新たに形成された筋細胞と筋線維の融合に必要である58,61。したがって、両患者群におけるこのプロセスの活性増加は、筋線維の不安定性によって引き起こされる損傷に対する修復反応を示唆している。
各ネマリンミオパチーの影響はよく相関しており（r = 0.736）、妥当な重複が見られた（補足図11A～B）。これは、ACTA1およびTNNT1ネマリンミオパチーがプロテオームに同様の影響を与えることを示している。しかし、一部のタンパク質はACTA1またはTNNT1ネマリンミオパチーでのみ調節されていた（補足図11AおよびC）。線維化促進タンパク質MFAP4はTNNT1ネマリンミオパチーで最も発現が上昇したタンパク質の1つであったが、ACTA1ネマリンミオパチーでは変化がなかった。HOX遺伝子転写の調節に関与するPAF1C複合体の構成要素であるSKIC8は、TNNT1ネマリンミオパチーでは発現が低下したが、ACTA1ネマリンミオパチーでは影響を受けなかった（補足図11A）。 ACTA1ネマリンミオパチーとTNNT1ネマリンミオパチーを直接比較したところ、TNNT1ネマリンミオパチーではミトコンドリアタンパク質の減少と免疫系タンパク質の増加がより顕著であることが明らかになった（図6G～Hおよび補足図11C、11H～I）。これらのデータは、TNNT1ネマリンミオパチーで観察された萎縮/ジストロフィーの程度がTNNT1ネマリンミオパチーよりも大きいことと一致しており（図6A）、TNNT1ネマリンミオパチーの方がより重篤な病型であることを示唆している。
ネマリンミオパチーの観察された影響が筋肉全体レベルで持続するかどうかを評価するために、同じコホートのTNNT1ネマリンミオパチー患者の筋生検のバルクプロテオミクス解析を行い、対照群（各群n=3）と比較した（補足図13A、補足データセット25）。予想どおり、対照群は主成分分析で密接に関連していたが、TNNT1ネマリンミオパチー患者は単一線維解析で見られたのと同様にサンプル間の変動が大きかった（補足図13B）。バルク解析は、個々の線維を比較することで強調された差次的発現タンパク質（補足図13C、補足データセット26）および生物学的プロセス（補足図13D、補足データセット27）を再現したが、異なる線維タイプを区別する能力を失い、線維間の不均一な疾患の影響を考慮できなかった。
これらのデータを総合すると、単一筋線維プロテオミクスは、イムノブロッティングなどの標的型解析法では検出できない臨床的生物学的特徴を解明できることが示される。さらに、これらのデータは、表現型適応を説明するためにアクチン線維型分類（MYH）のみを使用することの限界を浮き彫りにする。実際、線維型の切り替えはアクチンネマリンミオパチーとトロポニンネマリンミオパチーで異なるものの、どちらのネマリンミオパチーも、MYH線維型分類と骨格筋線維代謝を切り離し、より速く、より酸化の少ない筋プロテオームへと導く。
細胞の不均一性は、組織が多様な要求を満たすために不可欠です。骨格筋では、これはしばしば、異なる程度の力発生と疲労度によって特徴付けられる線維タイプとして説明されます。しかし、これは骨格筋線維の多様性のほんの一部しか説明しておらず、その多様性はこれまで考えられていたよりもはるかに多様で複雑かつ多面的であることは明らかです。技術の進歩により、骨格筋線維を制御する因子が明らかになってきました。実際、私たちのデータは、タイプ2X線維が明確な骨格筋線維サブタイプではない可能性を示唆しています。さらに、代謝タンパク質、リボソームタンパク質、および細胞関連タンパク質が骨格筋線維の不均一性の主要な決定因子であることを特定しました。線虫ミオパチー患者のサンプルに私たちのプロテオミクスワークフローを適用することで、MYHベースの線維型分類は、特にシステムが撹乱されている場合、骨格筋の不均一性を完全に反映していないことをさらに実証しました。実際、MYHを基盤とする線維の種類に関わらず、線虫性筋症は、より速く、酸化力の低い線維への変化をもたらす。
骨格筋線維は19世紀以来分類されてきました。最近のオミクス解析により、異なるMYH線維タイプの発現プロファイルと、さまざまな刺激に対するそれらの応答を理解し始めることができました。ここで説明するように、オミクスアプローチは、骨格筋線維タイプを定義するために単一（または少数）のマーカーの定量化に頼ることなく、従来の抗体ベースの方法よりも線維タイプマーカーを定量化する感度が高いという利点もあります。私たちは相補的なトランスクリプトミクスとプロテオミクスのワークフローを使用し、結果を統合して、ヒト骨格筋線維における線維の不均一性の転写および転写後調節を調べました。このワークフローの結果、健康な若い男性のコホートの外側広筋において、タンパク質レベルで純粋な2X型線維を特定できませんでした。これは、健康な外側広筋で純粋な2X線維が1%未満であることが判明した以前の単一線維研究と一致していますが、これは今後他の筋肉で確認する必要があります。 mRNA レベルではほぼ純粋な 2X 線維が検出され、タンパク質レベルでは混合 2A/2X 線維しか検出されないという矛盾は不可解である。MYH アイソフォーム mRNA の発現は概日リズムを示さないため、RNA レベルでは一見純粋な 2X 線維で MYH2 の開始シグナルを「見逃した」可能性は低い。考えられる説明の 1 つは、純粋に仮説ではあるが、MYH アイソフォーム間でタンパク質および/または mRNA の安定性に違いがあることである。実際、どの MYH アイソフォームに対しても 100% 純粋な速筋線維はなく、MYH1 mRNA の発現レベルが 70～90% の範囲にある場合、タンパク質レベルで MYH1 と MYH2 の量が等しくなるかどうかは不明である。しかし、トランスクリプトームまたはプロテオーム全体を考慮すると、クラスター分析では、正確な MYH 組成に関係なく、遅筋線維と速筋線維を表す 2 つの明確なクラスターのみを確実に識別できる。これは、通常2つの異なる筋核クラスターのみを同定する単一核トランスクリプトーム解析法を用いた分析結果と一致している。68, 69, 70 さらに、これまでのプロテオミクス研究では2X型線維が同定されているものの、これらの線維は他の速筋線維とは分離してクラスターを形成しておらず、MYHに基づく他の線維型と比較して、発現量が異なるタンパク質はごく少数しか示さない。14 これらの結果は、20世紀初頭の筋線維分類の考え方に戻るべきであることを示唆している。当時の分類では、ヒト骨格筋線維はMYHに基づく3つの異なるクラスではなく、代謝特性と収縮特性に基づく2つのクラスターに分類されていた。63
さらに重要なことに、筋線維の不均一性は複数の次元で考慮されるべきである。これまでの「オミクス」研究はこの方向性を示しており、骨格筋線維は離散的なクラスターを形成するのではなく、連続体に沿って配置されていることを示唆している。11, 13, 14, 64, 71 本研究では、骨格筋の収縮特性と代謝特性の違いに加えて、筋線維は細胞間相互作用と翻訳機構に関連する特徴によって区別できることを示す。実際、我々は骨格筋線維にリボソームの不均一性を発見し、これは遅筋線維と速筋線維のタイプとは無関係に不均一性に寄与している。遅筋線維と速筋線維の種類に関係なく、この大きな筋線維の不均一性の根本的な原因は依然として不明ですが、特定の力と負荷に最適に反応する筋束内の特殊な空間的組織、72 筋微小環境内の他の細胞型との特殊な細胞または器官特異的なコミュニケーション、73,74,75 または個々の筋線維内のリボソーム活性の違いを示している可能性があります。実際、RPL3 と RPL3L のパラログ置換または rRNA の 2'O-メチル化のレベルのいずれかによるリボソームヘテロプラスミーは、骨格筋肥大と関連していることが示されています 76,77。マルチオミクスおよび空間的アプリケーションと個々の筋線維の機能的特徴付けを組み合わせることで、マルチオミクスレベルでの筋生物学の理解がさらに進展します 78。
ネマリンミオパチー患者の単一筋線維のプロテオームを解析することで、骨格筋の臨床病態生理を解明するための単一筋線維プロテオミクスの有用性、有効性、および適用可能性も実証しました。さらに、ワークフローをグローバルプロテオミクス解析と比較することで、単一筋線維プロテオミクスがグローバル組織プロテオミクスと同じ深さの情報を提供し、線維間の不均一性と筋線維タイプを考慮することでこの深さを拡張できることを実証できました。ACTA1およびTNNT1ネマリンミオパチーで観察された線維タイプ比の予想される（ただし変動がある）違いに加えて、MYHを介した線維タイプの切り替えとは無関係に酸化および細胞外リモデリングも観察されました。 TNNT1ネマリンミオパチーでは線維化が以前から報告されている。19 しかし、我々の分析では、ACTA1およびTNNT1ネマリンミオパチー患者の筋線維において、筋形質膜修復機構に関与するアネキシンなどの細胞外分泌ストレス関連タンパク質のレベルが増加していることも明らかにし、この知見をさらに発展させた。57,58,59 結論として、ネマリンミオパチー患者の筋線維におけるアネキシンレベルの上昇は、重度に萎縮した筋線維を修復するための細胞応答を表している可能性がある。
本研究は、これまでに実施されたヒトの単一筋線維全筋オミクス解析としては最大規模ではあるものの、いくつかの限界がある。我々は、比較的小規模で均質な被験者群と単一の筋肉（外側広筋）から骨格筋線維を分離した。そのため、筋肉の種類や筋肉生理の極端な状態において、特定の線維集団が存在する可能性を排除することはできない。例えば、高度に訓練されたスプリンターや筋力系アスリート79、あるいは筋不活動期間66,80において、超高速線維（例えば、純粋な2X線維）のサブセットが出現する可能性を排除することはできない。さらに、被験者のサンプルサイズが限られていたため、線維タイプの比率は男女で異なることが知られているにもかかわらず、線維の多様性における性差を調査することができなかった。また、同一の筋線維または同一の被験者からのサンプルについて、トランスクリプトーム解析とプロテオーム解析を実施することもできなかった。我々や他の研究者が、オミクス解析を用いて単一細胞および単一筋線維の解析を最適化し、超低サンプル投入量を実現する（本稿でミトコンドリアミオパチー患者の筋線維の解析で実証したように）につれて、単一の筋線維内でマルチオミクス（および機能的）アプローチを組み合わせる機会が明らかになってくる。
総じて、我々のデータは骨格筋の多様性の転写および転写後調節因子を特定し、説明するものである。具体的には、古典的なMYHに基づく筋線維タイプの定義に関連する、骨格筋生理学における長年の定説に異議を唱えるデータを示す。我々は、この議論を再燃させ、最終的には骨格筋線維の分類と多様性に関する我々の理解を再考することを期待する。
14名の白人参加者（男性12名、女性2名）が、この研究への参加に自発的に同意した。本研究はゲント大学病院倫理委員会（BC-10237）の承認を受け、2013年ヘルシンキ宣言に準拠し、ClinicalTrials.gov（NCT05131555）に登録された。参加者の一般的な特徴は補足表1に示されている。口頭および書面によるインフォームドコンセントを得た後、参加者は研究への最終的な参加前に健康診断を受けた。参加者は若年（22～42歳）、健康（既往症なし、喫煙歴なし）、適度な身体活動レベルであった。最大酸素摂取量は、以前に記載されたように、ステップエルゴメーターを用いて身体能力を評価することで測定した。81
筋肉生検サンプルは、安静時および絶食状態で、14日間隔で3回採取されました。これらのサンプルはより大規模な研究の一部として採取されたため、参加者は生検の40分前にプラセボ（乳糖）、H1受容体拮抗薬（フェキソフェナジン540mg）、またはH2受容体拮抗薬（ファモチジン40mg）を摂取しました。これらのヒスタミン受容体拮抗薬は安静時の骨格筋のフィットネスに影響を与えないことを以前に実証しており81、品質管理プロットでは状態に関連したクラスタリングは観察されませんでした（補足図3および6）。実験日の前に48時間、標準化された食事（体重1kgあたり41.4kcal、炭水化物5.1g、タンパク質1.4g、脂肪1.6g）を摂取し、実験日の朝には標準化された朝食（炭水化物1.5g）を摂取した。局所麻酔（エピネフリンを含まない1%リドカイン0.5ml）下で、経皮的ベルグストローム吸引法を用いて大腿外側広筋から筋生検組織を採取した。82 筋組織サンプルは直ちにRNAlaterに包埋し、手動による筋線維解剖まで4℃で保存した（最長3日間）。
新鮮に分離した筋線維束を培養皿内の新鮮なRNAlater培地に移した。次に、実体顕微鏡と細いピンセットを使用して個々の筋線維を手動で解剖した。生検組織ごとに25本の線維を解剖し、生検組織の異なる領域から線維を選択することに特に注意を払った。解剖後、各線維を、不要なタンパク質とDNAを除去するためにプロテイナーゼKとDNase酵素を含む3μlの溶解バッファー（SingleShot Cell Lysis Kit、Bio-Rad）に穏やかに浸した。細胞溶解とタンパク質/DNAの除去は、短時間のボルテックス、マイクロ遠心分離機での液体の遠心分離、および室温でのインキュベーション（10分）によって開始した。その後、溶解液をサーマルサイクラー（T100、Bio-Rad）で37℃で5分間、75℃で5分間インキュベートし、その後、さらなる処理まで直ちに-80℃で保存した。
2 µl の筋線維ライセートから QuantSeq-Pool 3′ mRNA-Seq Library Prep Kit (Lexogen) を使用して Illumina 互換のポリアデニル化 RNA ライブラリーを調製した。詳細な方法は製造元のマニュアルに記載されている。プロセスは、逆転写による第一鎖 cDNA 合成から始まる。この過程で、サンプルのプールを確実にし、下流処理中の技術的変動を低減するために、ユニーク分子識別子 (UMI) とサンプル固有の i1 バーコードが導入される。次に、96 本の筋線維からの cDNA をプールし、磁気ビーズで精製した後、RNA を除去し、ランダムプライマーを使用して第二鎖合成を行う。ライブラリーを磁気ビーズで精製し、プール固有の i5/i7 タグを追加し、PCR 増幅する。最終精製ステップで Illumina 互換ライブラリーが生成される。各ライブラリープールの品質は、High Sensitivity Small Fragment DNA Analysis Kit (Agilent Technologies、DNF-477-0500) を使用して評価した。
Qubitによる定量に基づき、プールされたサンプルは等モル濃度（2 nM）でさらに混合された。得られたプールサンプルは、NovaSeq S2試薬キット（1 × 100ヌクレオチド）と2 nMのローディング濃度（4% PhiX）を用いて、NovaSeq 6000装置で標準モードでシーケンスされた。
当社のパイプラインは、Lexogen の QuantSeq Pool データ解析パイプライン (https://github.com/Lexogen-Tools/quantseqpool_analysis) に基づいています。データはまず、i7/i5 インデックスに基づいて bcl2fastq2 (v2.20.0) を使用してデマルチプレックスされました。次に、リード 2 は i1 サンプル バーコードに基づいて idemux (v0.1.6) を使用してデマルチプレックスされ、UMI 配列は umi_tools (v1.0.1) を使用して抽出されました。その後、リードは cutadapt (v3.4) を使用して複数回トリミングされ、短いリード (長さ 20 未満) またはアダプター配列のみで構成されるリードが除去されました。次に、リードは STAR (v2.6.0c) を使用してヒトゲノムにアライメントされ、BAM ファイルは SAMtools (v1.11) を使用してインデックス化されました。重複リードは umi_tools (v1.0.1) を使用して除去されました。最後に、アライメントカウントはSubread（v2.0.3）のfeatureCountsを使用して実行されました。品質管理は、パイプラインのいくつかの中間段階でFastQC（v0.11.9）を使用して実行されました。
以降のバイオインフォマティクス処理と可視化はすべて R (v4.2.3) で実行され、主に Seurat (v4.4.0) ワークフローが使用されました。83 そのため、個々の UMI 値とメタデータ マトリックスは Seurat オブジェクトに変換されました。全繊維の 30% 未満で発現している遺伝子は削除されました。低品質のサンプルは、最小閾値 1000 の UMI 値と 1000 個の検出された遺伝子に基づいて削除されました。最終的に、925 本の繊維がすべての品質管理フィルタリング手順を通過しました。UMI 値は Seurat SCTransform v2 メソッド 84 を使用して正規化され、検出されたすべての 7418 個の特徴が含まれ、参加者間の差は回帰によって除去されました。関連するすべてのメタデータは、補足データセット 28 にあります。