PyFeat 完整教學

Repository: https://github.com/RafsanjaniHub/PyFeat Stars: 97 | Fork: 32 | Language: Python | License: GPL-3.0 Tags: computational-biology, bioinformatics, genomics, proteomics Homepage: http://rafsanjani.pythonanywhere.com/ 最後更新: 2025-11-30

2. 核心架構

2.1 整體流程架構


flowchart TD
    A["FASTA 序列檔案
(DNA / RNA / Protein)"] --> B["read.py
讀取 FASTA + Label"]
    B --> C["generateFeatures.py
特徵擷取引擎"]

    C --> D1["pseudoKNC
k-mer 頻率"]
    C --> D2["k-Gap 特徵
monoMono ~ triDi
(8 種組合)"]
    C --> D3["Z-Curve
(x, y, z 軸)"]
    C --> D4["GC Content
AT/GC Ratio
Cumulative Skew"]

    D1 --> E["合併特徵矩陣
numpy array"]
    D2 --> E
    D3 --> E
    D4 --> E

    E --> F{"輸出選擇"}
    F -->|"fullDataset=1"| G1["fullDataset.csv
全部特徵"]
    F -->|"optimumDataset=1"| G2["selectedImportantFeatures.py
AdaBoost 特徵選擇"]
    F -->|"testDataset=1"| G3["testDataset.csv
獨立測試集"]

    G2 --> G4["optimumDataset.csv
篩選後特徵"]

    G1 --> H["runClassifiers.py
10 種 ML 分類器
k-Fold CV"]
    G4 --> H

    H --> I1["evaluationResults.txt"]
    H --> I2["auROC.png"]
    H --> I3["Accuracy_boxPlot.png"]

    G4 --> J["trainModel.py
訓練最終模型"]
    J --> K["dumpModel.pkl"]
    K --> L["evaluateModel.py
獨立驗證"]

2.2 特徵類型詳解

PyFeat 支援的特徵可分為四大類：

A. Pseudo k-tuple Nucleotide Composition (pseudoKNC; 偽 k 元組核苷酸組成)

計算序列中所有 1-mer、2-mer、3-mer 的出現次數。對於 DNA/RNA（4 字母表），kTuple=3 時產生 4 + 16 + 64 = 84 個特徵；對於蛋白質（20 字母表），產生 20 + 400 + 8000 = 8,420 個特徵。

B. k-Gap 特徵（8 種 mono/di/tri 組合）

以「左側子序列 + 間隔 + 右側子序列」的模式計數，例如 monoMonoKGap (X_X) 計算兩個單一殘基之間存在 k 個間隔的模式數。8 種組合為：

名稱模式DNA 特徵數 (kGap=1)蛋白質特徵數 (kGap=1)
monoMonoX_X16400
monoDiX_XX648,000
monoTriX_XXX256160,000
diMonoXX_X648,000
diDiXX_XX256160,000
diTriXX_XXX1,0243,200,000
triMonoXXX_X256160,000
triDiXXX_XX1,0243,200,000

C. 全域序列統計特徵（僅 DNA/RNA）

  • Z-Curve: x = (A+G) - (C+T), y = (A+C) - (G+T), z = (A+T) - (C+G)（3 個特徵）
  • GC Content: (G+C) / (A+C+G+T) x 100（1 個特徵）
  • Cumulative Skew: GC Skew = (G-C)/(G+C), AT Skew = (A-T)/(A+T)（2 個特徵）
  • AT/GC Ratio: (A+T) / (G+C)（1 個特徵）

D. 特徵選擇（AdaBoost-based）

使用 AdaBoost（n_estimators=500）的 feature_importances_ 屬性，保留重要性 > 0.00 的特徵，自動輸出 selectedIndex.txt 追蹤被選取的特徵索引。

2.3 內建分類器


flowchart LR
    subgraph classifiers["10 種內建分類器"]
        LR["Logistic Regression
(LR; 邏輯迴歸)"]
        KNN["k-Nearest Neighbor
(KNN; k 近鄰)"]
        DT["Decision Tree
(DT; 決策樹)"]
        NB["Naive Bayes
(NB; 樸素貝氏)"]
        BAG["Bagging
(袋裝法)"]
        RF["Random Forest
(RF; 隨機森林)"]
        AB["AdaBoost
(自適應增強)"]
        GB["Gradient Boosting
(GB; 梯度提升)"]
        SVM["Support Vector Machine
(SVM; 支持向量機)"]
        LDA["Linear Discriminant Analysis
(LDA; 線性判別分析)"]
    end

    CSV["optimumDataset.csv"] --> classifiers
    classifiers --> OUT1["evaluationResults.txt
Accuracy / auROC / auPR / F1 / MCC"]
    classifiers --> OUT2["auROC.png
ROC 曲線圖"]
    classifiers --> OUT3["Accuracy_boxPlot.png
準確率盒鬚圖"]

評估指標包含：Accuracy（準確率）、auROC、auPR（Average Precision; 平均精確率）、F1 Score、MCC（Matthews Correlation Coefficient; 馬修斯相關係數）、Sensitivity（靈敏度）、Specificity（特異度）、以及 Confusion Matrix（混淆矩陣）。

3. 安裝與設定

3.1 環境需求

項目版本需求
Python>= 3.5
NumPy>= 1.13.0（核心特徵生成）
scikit-learn>= 0.19.0（分類器 / 特徵選擇）
pandas>= 0.21.0（CSV 讀取）
matplotlib>= 2.1.0（ROC / BoxPlot 繪圖）

3.2 使用 uv 安裝（推薦）

 1# 建立專案目錄與虛擬環境
 2mkdir pyfeat-workspace && cd pyfeat-workspace
 3uv venv --python 3.11
 4source .venv/bin/activate
 5
 6# 安裝依賴
 7uv pip install numpy scikit-learn pandas matplotlib
 8
 9# 下載 PyFeat
10git clone https://github.com/RafsanjaniHub/PyFeat.git
11cd PyFeat/Codes

3.3 使用 conda 安裝

1conda create -n pyfeat python=3.11 numpy scikit-learn pandas matplotlib -y
2conda activate pyfeat
3
4git clone https://github.com/RafsanjaniHub/PyFeat.git
5cd PyFeat/Codes

3.4 驗證安裝

 1# 使用內建 demo 資料快速測試
 2python main.py \
 3  -seq=DNA \
 4  -full=1 \
 5  -fa=../Datasets/DNA/demoFASTA.txt \
 6  -la=../Datasets/DNA/demoLabel.txt \
 7  -ktuple=2 \
 8  -kgap=1 \
 9  -pseudo=1 \
10  -zcurve=1 \
11  -gc=1
12
13# 預期輸出：
14# Datasets fetching done.
15# Features extraction begins. Be patient! ...
16# Features extraction ends.
17# [Total extracted feature: XX]
18# Converting (full) CSV is begin.
19# Converting (full) CSV is end.

4. 使用方式與程式碼範例

4.1 範例一：DNA Promoter 序列特徵生成與分類

這個範例使用 PyFeat 內建的 DNA promoter 資料集，生成完整特徵集、執行特徵選擇，然後用 10 種分類器進行比較。

 1# ---- Step 1：生成特徵 ----
 2cd PyFeat/Codes
 3
 4python main.py \
 5  --sequenceType=DNA \
 6  --fullDataset=1 \
 7  --optimumDataset=1 \
 8  --fasta=../Datasets/DNA/FASTA.txt \
 9  --label=../Datasets/DNA/Label.txt \
10  --kTuple=3 \
11  --kGap=5 \
12  --pseudoKNC=1 \
13  --zCurve=1 \
14  --gcContent=1 \
15  --cumulativeSkew=1 \
16  --atgcRatio=1 \
17  --monoMono=1 \
18  --monoDi=1 \
19  --monoTri=1 \
20  --diMono=1 \
21  --diDi=1 \
22  --diTri=1 \
23  --triMono=1 \
24  --triDi=1
25
26# 產出：
27#   fullDataset.csv       — 全部特徵（數千列 x 數千行）
28#   optimumDataset.csv    — AdaBoost 篩選後的重要特徵
29#   selectedIndex.txt     — 被選取特徵的原始索引
30
31# ---- Step 2：分類器比較 ----
32python runClassifiers.py \
33  --nFCV=10 \
34  --dataset=optimumDataset.csv \
35  --auROC=1 \
36  --boxPlot=1
37
38# 產出：
39#   evaluationResults.txt — 10 種分類器的完整評估指標
40#   auROC.png            — ROC 曲線比較圖
41#   Accuracy_boxPlot.png — 準確率盒鬚圖

輸入檔案格式說明：

FASTA 檔案（FASTA.txt）格式：

1>promoter-1
2AGCCAGGCGAGATATGATCTATATCAATTTCTCATCTATAATGCTTTGTTAGTATCT
3CGTCGCCGACTTAATAAAGAGAGA
4>promoter-2
5CGGGCCTATAAGCCAGGCGAGATATGATCTATATCAATTTCTCATCTATAATGCTTT
6GTTAGTATCTCGTCGCCGACTTAA

Label 檔案（Label.txt）格式（每行一個類別標籤，與 FASTA 序列一一對應）：

11
21
30
40

4.2 範例二：蛋白質折疊分類（訓練 + 獨立測試）

使用蛋白質結構分類資料集（PDB1075 訓練 / PDB186 獨立測試），完成特徵生成、模型訓練、獨立驗證的完整流程。

 1# ---- Step 1：生成訓練集特徵（含特徵選擇） ----
 2python main.py \
 3  --sequenceType=Protein \
 4  --optimumDataset=1 \
 5  --fasta=../Datasets/Protein/PDB1075_FASTA.txt \
 6  --label=../Datasets/Protein/PDB1075_Label.txt \
 7  --kTuple=2 \
 8  --kGap=3 \
 9  --pseudoKNC=1 \
10  --monoMono=1 \
11  --monoDi=1 \
12  --diMono=1 \
13  --diDi=1
14
15# ---- Step 2：生成獨立測試集特徵 ----
16# 注意：必須使用相同的特徵參數
17python main.py \
18  --sequenceType=Protein \
19  --testDataset=1 \
20  --fasta=../Datasets/Protein/PDB186_independentFASTA.txt \
21  --label=../Datasets/Protein/PDB186_independentLabel.txt \
22  --kTuple=2 \
23  --kGap=3 \
24  --pseudoKNC=1 \
25  --monoMono=1 \
26  --monoDi=1 \
27  --diMono=1 \
28  --diDi=1
29
30# ---- Step 3：訓練模型 ----
31python trainModel.py \
32  --dataset=optimumDataset.csv \
33  --model=RF
34
35# 產出：dumpModel.pkl
36
37# ---- Step 4：獨立測試集驗證 ----
38python evaluateModel.py \
39  --optimumDatasetPath=optimumDataset.csv \
40  --testDatasetPath=testDataset.csv
41
42# 評估流程：
43#   1. 讀取 selectedIndex.txt 取得訓練時選定的特徵索引
44#   2. 用訓練集的 mean/std 對測試集做 Z-score 標準化
45#   3. 載入 dumpModel.pkl 進行預測
46#   4. 輸出 Accuracy / Sensitivity / Specificity / MCC 等指標

4.3 範例三：整合到 Python 腳本中的自訂流程

如果你想在自己的 Python 管線中直接使用 PyFeat 的特徵生成模組（不透過 CLI），可以直接 import 核心模組：

 1"""
 2自訂 PyFeat 特徵生成管線
 3適用場景：將 PyFeat 特徵整合至下游 ML/DL pipeline
 4"""
 5import sys
 6import numpy as np
 7from argparse import Namespace
 8
 9# 將 PyFeat/Codes 加入 Python path
10sys.path.insert(0, '/path/to/PyFeat/Codes')
11
12import read
13import generateFeatures
14import selectedImportantFeatures
15
16# ---- Step 1：讀取序列 ----
17fasta_path = '/path/to/PyFeat/Datasets/DNA/FASTA.txt'
18label_path = '/path/to/PyFeat/Datasets/DNA/Label.txt'
19
20X_seqs, Y_labels = read.fetchXY(fasta_path, label_path)
21print(f"載入 {len(X_seqs)} 條序列，{len(set(Y_labels))} 個類別")
22
23# ---- Step 2：組裝參數（模擬 argparse） ----
24args = Namespace(
25    sequenceType='DNA',
26    kTuple=3,
27    kGap=3,            # 降低 kGap 以減少特徵維度
28    pseudoKNC=1,
29    zCurve=1,
30    gcContent=1,
31    cumulativeSkew=1,
32    atgcRatio=1,
33    monoMono=1,
34    monoDi=1,
35    monoTri=0,          # 關閉高維特徵以節省記憶體
36    diMono=1,
37    diDi=0,
38    diTri=0,
39    triMono=0,
40    triDi=0,
41    fullDataset=0,
42    testDataset=0,
43    optimumDataset=0,
44)
45
46# ---- Step 3：生成特徵 ----
47T = generateFeatures.gF(args, X_seqs, Y_labels)
48X_features = T[:, :-1]   # 特徵矩陣
49Y_target = T[:, -1]       # 標籤向量
50
51print(f"特徵矩陣維度：{X_features.shape}")
52# 範例輸出：特徵矩陣維度：(244, 492)
53
54# ---- Step 4：AdaBoost 特徵選擇 ----
55X_selected = selectedImportantFeatures.importantFeatures(
56    X_features, Y_target
57)
58print(f"篩選後特徵數：{X_selected.shape[1]}")
59
60# ---- Step 5：接續你的下游管線 ----
61# 例如：送入 XGBoost / LightGBM / 深度學習模型
62from sklearn.model_selection import train_test_split
63from sklearn.ensemble import GradientBoostingClassifier
64from sklearn.metrics import classification_report
65
66X_train, X_test, y_train, y_test = train_test_split(
67    X_selected, Y_target, test_size=0.2, stratify=Y_target, random_state=42
68)
69
70clf = GradientBoostingClassifier(n_estimators=200, max_depth=5)
71clf.fit(X_train, y_train)
72y_pred = clf.predict(X_test)
73
74print(classification_report(y_test, y_pred, digits=4))

5. 在蛋白質設計/基因組模擬生態系中的定位

5.1 生態系定位

在 Domain 5（蛋白質設計與基因組模擬）的 22 個專案中，PyFeat 扮演的是**基礎特徵工程層（Feature Engineering Layer）**的角色。它不直接做蛋白質設計或基因組模擬，而是將生物序列轉換為下游 ML 任務可用的數值表徵（Numerical Representations）。


flowchart TB
    subgraph upstream["上游：序列來源"]
        A1["基因組模擬
wgsim / HyenaDNA"]
        A2["蛋白質設計
RFdiffusion / Chroma / EvoDiff"]
        A3["單細胞模擬
scGPT / scvi-tools / splatter"]
        A4["實驗序列
NGS / Mass Spec"]
    end

    subgraph pyfeat["PyFeat：特徵工程層"]
        B1["DNA/RNA 特徵
pseudoKNC + k-Gap + Z-Curve"]
        B2["蛋白質特徵
pseudoKNC + k-Gap
(20 字母表)"]
    end

    subgraph downstream["下游：ML 應用"]
        C1["序列分類
Promoter / Enhancer 預測"]
        C2["蛋白質功能預測
折疊分類 / 結合位點"]
        C3["藥物靶點篩選
Target Identification"]
        C4["合成序列品質評估
Generated vs Real"]
    end

    A1 --> B1
    A2 --> B2
    A3 --> B1
    A4 --> B1
    A4 --> B2

    B1 --> C1
    B1 --> C3
    B2 --> C2
    B2 --> C3
    B1 --> C4
    B2 --> C4

5.2 與 AIKT Pipeline 的對應

AIKT WP 階段PyFeat 應用場景
WP1 靶點發現 (Target Discovery)將候選基因序列轉為特徵向量，用 ML 預測是否為有效靶點
WP2 靶點驗證 (Target Validation)比較 scGPT/scvi-tools 模擬的單細胞資料中，靶點相關基因的序列特徵分佈
WP4 先導物優化 (Lead Optimization)對蛋白質藥物（如抗體序列）做特徵化，輔助 binding affinity 預測
WP7 生物製劑 (Biologics)將 RFdiffusion/EvoDiff 生成的蛋白質序列特徵化，評估合成序列品質

5.3 藍海機會

  1. 合成序列品質評估器 (Synthetic Sequence Quality Assessor)：用 PyFeat 的特徵空間比較生成式模型（EvoDiff、Chroma）產出的序列與天然序列的分佈差異，建立品質指標（QC Score）
  2. 跨域特徵融合 (Cross-Domain Feature Fusion)：結合 PyFeat 的序列組成特徵與 ESM/ProstT5 的深度學習嵌入（Deep Learning Embeddings），建立 hybrid representation 提升下游預測效能
  3. 抗體序列優化回饋迴路 (Antibody Optimization Feedback Loop)：在 DiffAb 生成 → PyFeat 特徵化 → ML 篩選 → DiffAb 微調的迴路中，PyFeat 可作為快速篩選閘門

6. 與其他工具的整合

6.1 與 ESM（Evolutionary Scale Modeling）互補

ESM 產出的是高維密集嵌入（Dense Embeddings; 640-5120 維），PyFeat 產出的是稀疏可解釋特徵（Sparse Interpretable Features）。兩者可以互補：

 1"""
 2ESM 嵌入 + PyFeat 特徵的 hybrid representation
 3"""
 4import numpy as np
 5from sklearn.preprocessing import StandardScaler
 6
 7# 假設已有 ESM 嵌入和 PyFeat 特徵
 8esm_embeddings = np.load('esm_embeddings.npy')      # shape: (N, 1280)
 9pyfeat_features = np.load('pyfeat_features.npy')     # shape: (N, 500)
10
11# 標準化後拼接
12scaler_esm = StandardScaler()
13scaler_pf = StandardScaler()
14
15X_esm = scaler_esm.fit_transform(esm_embeddings)
16X_pf = scaler_pf.fit_transform(pyfeat_features)
17
18X_hybrid = np.hstack([X_esm, X_pf])
19print(f"Hybrid 特徵維度：{X_hybrid.shape}")
20# Hybrid 特徵維度：(N, 1780)

6.2 與 scGPT/scvi-tools 搭配

在單細胞分析中，可將差異表現基因（DEGs; Differentially Expressed Genes）的核苷酸序列匯入 PyFeat，產生序列層面的特徵補充 expression-level 的分析：

步驟工具輸入輸出
1. 單細胞資料分析scGPT / scvi-toolsscRNA-seq count matrixDEG 清單
2. 取得基因序列Ensembl API / BioMartGene IDsFASTA 序列
3. 特徵生成PyFeatFASTA + Labels數值特徵矩陣
4. 分類/預測scikit-learn / XGBoost特徵矩陣靶點優先排序

6.3 與 RFdiffusion/EvoDiff 搭配

用於評估生成式蛋白質設計工具的輸出品質：

步驟工具說明
1. 生成蛋白質RFdiffusion / EvoDiff產出候選蛋白質序列
2. 天然序列參考UniProt / PDB取得同功能家族的天然序列
3. 特徵擷取PyFeat對生成序列與天然序列分別生成特徵
4. 分佈比較scipy / scikit-learnKS test / t-SNE / UMAP 視覺化特徵空間差異

6.4 與 paper-qa-lite / paper-search 搭配（AIKT 管線內）

1paper-search: "sequence feature extraction bioinformatics k-mer k-gap"
2  → 找到相關方法論文獻
3  → PyFeat 實作特徵擷取
4  → paper-qa-lite: 針對產出的特徵矩陣問「哪些特徵對 X 預測最重要？」

7. 優缺點分析

7.1 優點

面向說明
統一介面 (Unified Interface)單一 CLI 同時處理 DNA / RNA / Protein，不需為不同序列類型切換工具
零依賴特徵生成核心特徵生成只需 NumPy，不需額外的生物資訊套件
端到端管線 (End-to-End Pipeline)從 FASTA → 特徵 → 特徵選擇 → 分類 → 評估一站完成
可解釋特徵 (Interpretable Features)k-mer / k-gap 特徵具有明確的生物學意義，可直接追溯到序列模式
內建分類器比較10 種分類器 + 自動 ROC / BoxPlot，適合快速基準測試（Benchmarking）
輕量級 (Lightweight)無 GPU 需求、無深度學習依賴，適合計算資源有限的環境

7.2 缺點

面向說明建議替代/補充
無深度學習表徵僅提供傳統統計特徵，無法捕捉長程序列依賴性（Long-range Dependencies）搭配 ESM / ProstT5 / HyenaDNA 的嵌入層
蛋白質高維爆炸蛋白質 k-gap 特徵可達百萬維，記憶體消耗極大降低 kGap / kTuple 參數，或僅啟用部分特徵類型
無 pip 安裝未打包為 PyPI 套件，需手動 clone 並調整 sys.path自行包裝或固定版本的 git submodule
檔案 I/O 耦合輸出直接寫入工作目錄的固定檔名（如 fullDataset.csv），不支援自訂路徑修改 save.py 或在外層腳本中管理工作目錄
無增量處理每次執行都從頭計算所有特徵，無快取（Cache）機制可自行加入 joblib.Memory 快取層
維護頻率低最後更新為 2025-11-30，無持續性的版本發佈功能穩定但不會有新特徵加入
僅支援二分類內建分類器和評估指標預設為二元分類（Binary Classification）多分類任務需自行修改 runClassifiers.py

7.3 適用場景總結

適合不適合
快速基準測試：比較不同序列特徵對分類效果的影響需要深度學習端到端訓練的場景
可解釋性需求高的生物資訊分析超長序列（>10K bp）的全基因組分析
計算資源有限、不需 GPU 的環境需要即時更新或串流處理的管線
蛋白質/DNA 分類任務的特徵工程前處理大規模蛋白質組學（Proteomics）高通量分析
教學用途：展示序列特徵擷取的原理生產環境需要 API 化部署的場景

附錄：快速參考卡

 1# ========== 特徵生成 ==========
 2python main.py -seq=DNA -full=1 -optimum=1 \
 3  -fa=FASTA.txt -la=Label.txt \
 4  -ktuple=3 -kgap=5 \
 5  -pseudo=1 -zcurve=1 -gc=1 -skew=1 -atgc=1 \
 6  -f11=1 -f12=1 -f13=1 -f21=1 -f22=1 -f23=1 -f31=1 -f32=1
 7
 8# ========== 分類器比較 ==========
 9python runClassifiers.py -cv=10 -data=optimumDataset.csv -roc=1 -box=1
10
11# ========== 訓練模型 ==========
12python trainModel.py -data=optimumDataset.csv -m=RF
13
14# ========== 獨立測試集驗證 ==========
15python evaluateModel.py -optimumPath=optimumDataset.csv -testPath=testDataset.csv
16
17# ========== 可用分類器代碼 ==========
18# LR | SVM | KNN | DT | NB | Bagging | RF | AB | GB | LDA

文件產生時間: 2026-06-10 | 資料來源: GitHub RafsanjaniHub/PyFeat (commit HEAD) | 適用版本: PyFeat master branch