scikit-bio 是做什麼的?
概述
scikit-bio 是一個全面的 Python 函式庫,專為處理生物資料而設計。在需要進行生物資訊學分析時,包括序列操作、比對、系統發育學、微生物生態學和多元統計,請應用此技能。
何時使用此技能
當使用者符合以下情況時,應使用此技能:
- 處理生物序列 (DNA、RNA、蛋白質)
- 需要讀取/寫入生物檔案格式 (FASTA、FASTQ、GenBank、Newick、BIOM 等)
- 執行序列比對或搜尋 motif
- 建構或分析系統發育樹
- 計算多樣性指標 (alpha/beta 多樣性、UniFrac 距離)
- 執行排序分析 (PCoA、CCA、RDA)
- 對生物/生態資料執行統計檢定 (PERMANOVA、ANOSIM、Mantel)
- 分析微生物組或群落生態學資料
- 處理來自語言模型的蛋白質嵌入
- 需要操作生物資料表格
核心功能
1. 序列操作
使用專門的 DNA、RNA 和蛋白質資料類別來處理生物序列。
主要操作:
- 從 FASTA、FASTQ、GenBank、EMBL 格式讀取/寫入序列
- 序列切片、連接和搜尋
- 反向互補、轉錄 (DNA→RNA) 和轉譯 (RNA→蛋白質)
- 使用正則表達式尋找 motif 和模式
- 計算距離 (Hamming、k-mer 爲基礎)
- 處理序列品質分數和中繼資料
常見模式:
import skbio
# 從檔案讀取序列
seq = skbio.DNA.read('input.fasta')
# 序列操作
rc = seq.reverse_complement()
rna = seq.transcribe()
protein = rna.translate()
# 尋找 motif
motif_positions = seq.find_with_regex('ATG[ACGT]{3}')
# 檢查屬性
has_degens = seq.has_degenerates()
seq_no_gaps = seq.degap()
重要注意事項:
- 對於需要驗證的語法序列,使用
DNA、RNA、Protein類別 - 對於沒有字母限制的通用序列,使用
Sequence類別 - 品質分數會自動從 FASTQ 檔案載入到位元位置中繼資料中
- 中繼資料類型:序列層級 (ID、描述)、位置層級 (每個鹼基)、區間層級 (區域/特徵)
2. 序列比對
使用 pair_align 引擎 (於 scikit-bio 0.7.0 引入) 執行成對和多序列比對,這是一個多功能且高效的動態規劃比對器。
主要功能:
- 在一個函數中進行全域、局部和半全域比對 (自由端點可配置)
- 便捷包裝函數
pair_align_nucl(類似 BLASTN) 和pair_align_prot(類似 BLASTP) - 可配置的評分:匹配/不匹配元組或命名的取代矩陣;線性或仿射空位罰分
PairAlignPath結果攜帶 CIGAR 字串,並可轉換爲比對後的序列- 使用
TabularMSA儲存和操作多序列比對
常見模式:
from skbio import DNA, Protein
from skbio.alignment import pair_align_nucl, pair_align_prot, pair_align, TabularMSA
# 使用 BLASTN 預設值的核酸比對
seq1, seq2 = DNA('ACTACCAGATTACTTACGGATCAGG'), DNA('CGAAACTACTAGATTACGGATCTTA')
aln = pair_align_nucl(seq1, seq2)
aln.score # 比對分數 (浮點數)
path = aln.paths[0] # PairAlignPath (repr 顯示 CIGAR)
aligned_seqs = path.to_aligned((seq1, seq2)) # 有間隔字串的清單
# 從比對路徑和原始序列建構 TabularMSA
msa = TabularMSA.from_path_seqs(path, (seq1, seq2))
# 透過 pair_align 自訂演算法 (預設 mode='global')
aln = pair_align(seq1, seq2, mode='local') # Smith-Waterman
aln = pair_align(seq1, seq2, sub_score=(2, -3), gap_cost=(5, 2)) # 仿射空位
aln = pair_align(seq1, seq2, sub_score='NUC.4.4', gap_cost=3) # 取代矩陣,線性空位
# 蛋白質比對 (類似 BLASTP,BLOSUM62)
aln = pair_align_prot(Protein('HEAGAWGHEE'), Protein('PAWHEAE'))
# 從檔案讀取多序列比對並進行摘要
msa = TabularMSA.read('alignment.fasta', constructor=DNA)
consensus = msa.consensus()
重要注意事項:
pair_align取代了已移除的 SSW 包裝函數 (local_pairwise_align_ssw、StripedSmithWaterman) 和已棄用的純 Python 比對器 (global_pairwise_align、local_pairwise_align_nucleotide等)- 結果是
PairAlignResult,也可以解包爲score, paths, matrices(使用keep_matrices=True保留 DP 矩陣) sub_score接受(match, mismatch)元組或矩陣名稱 (例如'NUC.4.4'、'BLOSUM62');gap_cost接受單個數字 (線性) 或(open, extend)元組 (仿射)- 使用
PairAlignPath.from_cigar('1I8M2D5M2I')解析外部的 CIGAR 字串;使用align_score(...)對現有比對進行評分,使用align_dists(...)從 MSA 建構距離矩陣
3. 系統發育樹
建構、操作和分析代表演化關係的系統發育樹。
主要功能:
- 從距離矩陣建構樹 (UPGMA/WPGMA、Neighbor Joining、GME、BME)
- 使用最近鄰交換 (
nni) 進行樹重排 - 樹操作 (修剪、重新定根、遍歷)
- 距離計算 (透過
cophenet計算 patristic 距離,透過compare_rfd計算 Robinson-Foulds 距離) - ASCII 視覺化
- Newick 格式 I/O
常見模式:
from skbio import TreeNode
from skbio.tree import nj, upgma, gme, bme, rf_dists
# 從檔案讀取樹
tree = TreeNode.read('tree.nwk')
# 從距離矩陣建構樹
tree = nj(distance_matrix)
# 樹操作
subtree = tree.shear(['taxon1', 'taxon2', 'taxon3'])
tips = [node for node in tree.tips()]
lca = tree.lca(['taxon1', 'taxon2'])
# 計算距離
patristic_dist = tree.find('taxon1').distance(tree.find('taxon2'))
cophenetic_dm = tree.cophenet() # 末端節點間的 patristic 距離矩陣
# 比較兩棵樹 (Robinson-Foulds)
rf_distance = tree.compare_rfd(other_tree)
# 多棵樹間的成對 RF 距離 -> DistanceMatrix
rf_dm = rf_dists([tree, other_tree, third_tree])
重要注意事項:
- 使用
nj()進行鄰接法 (經典的系統發育方法) - 使用
upgma()進行 UPGMA/WPGMA (假設分子時鐘) - GME 和 BME 對於大型樹具有高度可擴展性;使用
nni()優化拓撲 cophenet()(原名tip_tip_distances) 返回 patristic 距離矩陣;compare_rfd()是 Robinson-Foulds 方法 (compare_wrfd/compare_cophenet適用於加權/cophenetic 變體)lca()是最低共同祖先;lowest_common_ancestor保留爲別名- 樹可以有根或無根;某些指標需要特定的定根方式
4. 多樣性分析
計算微生物生態學和群落分析的 alpha 和 beta 多樣性指標。
主要功能:
- Alpha 多樣性:豐富度 (
sobs、observed_features、chao1、ace)、Shannon、Simpson、Hill 數 (hill)、Faith's PD (faith_pd)、廣義 PD (phydiv)、Pielou's evenness - Beta 多樣性:Bray-Curtis、Jaccard、加權/未加權 UniFrac、歐幾里得距離
- 系統發育多樣性指標 (需要輸入樹)
- 稀疏化和子抽樣
- 與排序和統計檢定整合
常見模式:
from skbio.diversity import alpha_diversity, beta_diversity
# Alpha 多樣性 (系統發育指標需要 taxa= 來進行末端節點名稱映射)
alpha = alpha_diversity('shannon', counts_matrix, ids=sample_ids)
faith_pd = alpha_diversity('faith_pd', counts_matrix, ids=sample_ids,
tree=tree, taxa=feature_ids)
# Beta 多樣性
bc_dm = beta_diversity('braycurtis', counts_matrix, ids=sample_ids)
unifrac_dm = beta_diversity('unweighted_unifrac', counts_matrix,
ids=sample_ids, tree=tree, taxa=feature_ids)
# 取得可用的指標
from skbio.diversity import get_alpha_diversity_metrics
print(get_alpha_diversity_metrics())
重要注意事項:
- 計數必須是代表豐度的整數,而非相對頻率
- 系統發育指標的參數是
taxa=(在 0.6.0 中從otu_ids重新命名;舊名稱是已棄用的別名);observed_otus現在是observed_features(或sobs) counts_matrix可以是任何類似表格的輸入 (NumPy 陣列、pandas/polars DataFrame、BIOMTable或 AnnData),透過分派系統處理- 系統發育指標 (Faith's PD、UniFrac) 需要樹和分類單元到末端節點的映射
- 對於特定的樣本配對,使用
partial_beta_diversity();對於大型區塊分解計算,使用block_beta_diversity() - Alpha 多樣性返回
pandas.Series,beta 多樣性返回DistanceMatrix
5. 排序方法
將高維度生物資料降維到可視化的低維度空間。
主要功能:
- PCoA (主座標分析) 從距離矩陣
- CA (對應分析) 用於列聯表
- CCA (典型對應分析) 帶有環境約束
- RDA (冗餘分析) 用於線性關係
- 雙標圖投影用於特徵解釋
常見模式:
from skbio.stats.ordination import pcoa, cca
import skbio
# 從距離矩陣進行 PCoA (大型矩陣限制維度)
pcoa_results = pcoa(distance_matrix, dimensions=3)
pc1 = pcoa_results.samples['PC1']
pc2 = pcoa_results.samples['PC2']
# 內建的散點圖,按中繼資料欄位著色
fig = pcoa_results.plot(sample_metadata, column='bodysite')
# 帶有環境變數的 CCA
cca_results = cca(species_matrix, environmental_matrix)
# 儲存/載入排序結果
pcoa_results.write('ordination.txt')
results = skbio.OrdinationResults.read('ordination.txt')
重要注意事項:
- PCoA 適用於任何距離/相異性矩陣;將
dimensions作爲整數 (計數) 或 (0, 1] 範圍內的浮點數 (保留的累積變異比例) 傳遞 OrdinationResults暴露基於 pandas 的屬性:samples、features、eigvals、proportion_explained、biplot_scores、sample_constraints- CCA 揭示環境對群落組成的驅動因素
OrdinationResults.plot()產生 matplotlib 圖形;結果也可以與 seaborn/plotly 整合
6. 統計檢定
執行針對生態和生物資料的假設檢定。
主要功能:
- PERMANOVA:使用距離矩陣檢驗組間差異
- ANOSIM:組間差異的替代檢定
- PERMDISP:檢驗組間離散度的同質性
- Mantel 檢定:距離矩陣間的相關性
- Bioenv:找尋與距離相關的環境變數
- 差異豐度:
skbio.stats.composition中的ancom、dirmult_ttest和dirmult_lme(縱向混合效應)
常見模式:
from skbio.stats.distance import permanova, anosim, mantel
# 檢驗組間是否存在顯著差異
permanova_results = permanova(distance_matrix, grouping, permutations=999)
print(f"p-value: {permanova_results['p-value']}")
# ANOSIM 檢定
anosim_results = anosim(distance_matrix, grouping, permutations=999)
# 兩個距離矩陣間的 Mantel 檢定
mantel_results = mantel(dm1, dm2, method='pearson', permutations=999)
print(f"Correlation: {mantel_results[0]}, p-value: {mantel_results[1]}")
# 對特徵表格進行差異豐度分析 (建議使用原始計數)
from skbio.stats.composition import dirmult_ttest
da = dirmult_ttest(counts_table, grouping, treatment='caseA', reference='control')
重要注意事項:
- 置換檢定提供了非參數顯著性檢驗
- 使用 999 次或更多次置換以獲得穩健的 p 值
- PERMANOVA 對離散度差異敏感;建議與 PERMDISP 搭配使用
- Mantel 檢定評估矩陣相關性 (例如,地理距離與遺傳距離)
- 對差異豐度檢定提供原始計數,而非預先正規化的比例,以保留量級資訊
7. 檔案 I/O 與格式轉換
讀取和寫入 19 種以上生物檔案格式,並具備自動格式偵測功能。
支援格式:
- 序列:FASTA、FASTQ、GenBank、EMBL、QSeq
- 比對:Clustal、PHYLIP、Stockholm
- 樹:Newick
- 表格:BIOM (HDF5 和 JSON)
- 距離:帶分隔符號的方陣
- 分析:BLAST+6/7、GFF3、排序結果
- 中繼資料:經過驗證的 TSV/CSV
常見模式:
import skbio
# 使用自動格式偵測讀取
seq = skbio.DNA.read('file.fasta', format='fasta')
tree = skbio.TreeNode.read('tree.nwk')
# 寫入檔案
seq.write('output.fasta', format='fasta')
# 大型檔案的產生器 (節省記憶體)
for seq in skbio.io.read('large.fasta', format='fasta', constructor=skbio.DNA):
process(seq)
# 轉換格式
seqs = list(skbio.io.read('input.fastq', format='fastq', constructor=skbio.DNA))
skbio.io.write(seqs, format='fasta', into='output.fasta')
重要注意事項:
- 對於大型檔案,使用產生器以避免記憶體問題
- 當指定
into參數時,格式可以自動偵測 - 某些物件可以寫入多種格式
- 支援使用
verify=False進行標準輸入/輸出管道
8. 距離矩陣
使用統計方法建立和操作距離/相異性矩陣。
主要功能:
- 儲存對稱 (
DistanceMatrix,對角線爲零) 或一般成對 (PairwiseMatrix) 資料 - 基於 ID 的索引和切片
- 與多樣性、排序和統計檢定整合
- 讀取/寫入帶分隔符號的文字格式
常見模式:
from skbio import DistanceMatrix
import numpy as np
# 從陣列建立
data = np.array([[0, 1, 2], [1, 0, 3], [2, 3, 0]])
dm = DistanceMatrix(data, ids=['A', 'B', 'C'])
# 存取距離
dist_ab = dm['A', 'B']
row_a = dm['A']
# 從檔案讀取
dm = DistanceMatrix.read('distances.txt')
# 在下游分析中使用
pcoa_results = pcoa(dm)
permanova_results = permanova(dm, grouping)
重要注意事項:
DistanceMatrix強制對稱性和零 (空心) 對角線;它是SymmetricMatrix的子類別PairwiseMatrix(從DissimilarityMatrix重新命名,後者保留爲已棄用的別名) 允許一般/非對稱值- ID 實現與中繼資料和生物知識的整合