Lecture Notes: Introduction to Bioinformatics

AI3073 Introduction to Bioinformatics

Dr. Jiaxing CHEN

Chapter 2 Part 1: Gene and RNA

1. 基本概念 (Fundamentals)

• Gene (基因): A locus/region of DNA encoding a functional protein or RNA product; the molecular unit of heredity (遗传单位)。

• Central Dogma (中心法则): DNA → RNA → Protein。

  • Transcription (转录): DNA → pre-mRNA
  • Splicing (剪接): pre-mRNA → mRNA (去掉 introns, 保留 exons)
  • Translation (翻译): mRNA → Protein

---

2. Genome & 基因组特征

• Human genome length: ~3.2 billion base pairs (bp)。
• Cell scale: 从病毒 (nm) 到 human cell (~10 μm)，范围极大 (见 page 10 diagram)。
• Codon Table (密码子表, page 12): 3 个碱基 = 1 个 amino acid. Stop codons: UAA, UAG, UGA。

---

3. Computational Gene Prediction (计算基因预测)

• Challenge (难点):

  • Mammalian genomes 很大 (~3 billion bp)
  • <2% DNA 编码蛋白
  • Non-coding RNAs 难预测。

• Approaches (方法):

  1. Ab initio: 根据统计特征 (如 codon usage, ORF)
  2. Homology-based: 基于已知基因相似性
  3. Hybrid: 混合方法

• Key features (特征):

  • ORF (开放阅读框): start codon ATG → stop codon
  • Codon Usage 偏好
  • Motifs (启动子, enhancers, UTRs)

---

4. RNA 分类 (RNA Types)

• mRNA: 编码蛋白

• Non-coding RNA（ncRNA, 非编码RNA）：

  • tRNA: adaptor, 携带氨基酸 (~80 nt, clover-leaf structure)
  • rRNA: 组成 ribosome，催化蛋白合成
  • miRNA (~22 nt): 基因表达调控 (silencing), 有发夹状 precursor
  • lncRNA (>200 nt): 调控作用, scaffold/guide/enhancer
  • circRNA: 环状，稳定，不易降解，可来源于蛋白编码基因

---

5. RNA Secondary Structure (RNA 二级结构)

• Levels (层次结构):

  • Primary structure (一级结构): 核苷酸序列
  • Secondary structure (二级结构): 发夹、茎环等 (page 49)
  • Tertiary structure (三级结构): 二级结构之间的相互作用 (page 47)

• Prediction (预测方法):

  1. Co-variation analysis: 多序列比对, conserved base pairs
  2. Single sequence prediction: 最小自由能 (minimum free energy, MFE)

---

6. 核心 takeaway

• Genes = DNA functional unit
• Central Dogma: DNA → RNA → Protein
• Gene prediction: 难，因为非编码区很多
• ncRNAs: 种类繁多 (tRNA, rRNA, miRNA, lncRNA, circRNA)，功能复杂
• RNA structure: 从一级到三级，预测依赖计算方法

---

Chapter 2 Part 2: Genome, NGS, Transcriptome and Assembly

1. Overview

本章主要涵盖四个主题：

1. Genome (基因组)
2. DNA Sequencing (DNA 测序)
3. Transcriptome (转录组)
4. Genome Assembly (基因组组装)

---

2. Genome 基因组

• Definition (定义):
  A genome is the complete set of DNA (包括所有基因和非编码区) in an organism.
  → 它包含所有遗传信息，用于指导细胞功能和发育。

• Genetic Variation (遗传变异):

  • SNP (Single Nucleotide Polymorphism): 单个碱基变化
  • Indel (Insertion/Deletion): 小片段插入或缺失
  • CNV (Copy Number Variation): 拷贝数变化
  • SV (Structural Variation): 大规模结构变异（如倒位、易位等）
    ⭐ Calling genetic variation = 检测这些变异的过程。

---

3. DNA Sequencing (DNA测序)

From Sanger to NGS

• Sanger Sequencing: 第一代测序，基于链终止法（Dideoxy method）。

• Next Generation Sequencing (NGS, 下一代测序):

  • 高通量（High-throughput）
  • 低成本（Low cost per base）
  • 可同时读取数百万条 DNA 片段。

Sequencing Output Formats

• FASTA format:

  >sequence_name
  ATCGTAGCTA...

  ➜ 只包含序列。

• FASTQ format:
  包含序列 + 每个碱基的质量分值 (Phred score)。

  @sequence_name
  ATCGTAGCTA
  +
  FJJJJJJGJH

  ➜ 常用于 NGS 原始数据。

---

4. Transcriptome 转录组

• Definition (定义):
  所有转录出来的 RNA 分子的总和 (包括 mRNA + 非编码RNA)。

• RNA-seq (转录组测序):

  • 定量测量基因表达水平 (quantify transcript abundance)
  • 可发现新转录本和可变剪接事件。

Special Types

• Single-cell RNA sequencing (scRNA-seq)
  单细胞水平研究基因表达差异，揭示细胞异质性 (heterogeneity)。
• Spatial Transcriptomics (空间转录组学)
  同时保留 RNA 表达与空间位置信息，理解组织结构中的基因表达格局。

---

5. Genome Assembly 基因组组装

Purpose

将短测序片段 (short reads) 拼接成完整的基因组序列。

Types of Assembly

1. De novo assembly (从头组装):
   无参考基因组的情况下，通过片段重叠关系拼接。
2. Reference-based assembly (参考组装):
   将 reads 对齐到已有参考基因组。

Steps

1. Pre-processing: 过滤低质量 reads
2. Overlap / de Bruijn Graph Construction: 构建序列关系
3. Contig Assembly: 拼接连续序列片段 (contigs)
4. Scaffolding: 将 contigs 排列到染色体水平
5. Error Correction / Polishing: 提高准确度

---

6. Summary 核心总结

Topic	Key Concept	中文要点
Genome	Complete DNA content	生物的全部遗传信息
Variation Calling	SNP, CNV, SV	发现基因变异
NGS	High-throughput sequencing	高通量、低成本测序
Transcriptome	All RNA transcripts	表达调控研究核心
Assembly	Reconstruct genome	从 reads 拼出完整基因组

---

Chapter 2 Part 3: Alignment and NGS Reads Mapping - BLAST

Part 1. Biological Motivation & Sequence Alignment Fundamentals

（生物动机与序列比对基础）

1. The Biological Question

Question: How can we determine the similarity between two sequences?
我们想知道两个生物序列（如 DNA、RNA 或蛋白质）是否相似。

为什么重要？

• 相似序列 → 相似结构 → 相似功能 （Sequence → Structure → Function Paradigm）
• 相似序列 → 共同祖先（Common ancestor, Homology）

所以，“序列相似性”是研究基因功能与进化关系的出发点。

---

2. What is Sequence Alignment?

Definition: Sequence alignment 是一种计算方法，用于排列序列中对应的字符（碱基或氨基酸），以最大化它们的相似性。

目的：

将输入序列（inputted sequences）中的所有残基对齐，使得比对结果在功能或进化关系上达到最大相似度。

如 slide 中的例子（page 5），我们对齐两个序列 A 和 B，使得保守区域被正确地对应起来。

---

3. Types of Sequence Alignment

在生物信息学中常见两种主要类型的比对：

类型	名称	特点	应用场景
Global Alignment	全局比对	比对整个序列（end-to-end）	序列长度相似，整体比较
Local Alignment	局部比对	只比对局部最相似的部分	序列差异较大，仅部分保守

---

4. The Scoring System 评分系统

为了评估一个比对的好坏，需要定义一个 scoring function（评分函数）。
它由两部分组成：

1. Substitution Score（替换得分）
   当两个字符相同 → 得分高；不同 → 扣分。
   在蛋白质比对中，常用：

   • PAM 矩阵
   • BLOSUM 矩阵

   对 DNA 序列，可使用如下简单矩阵（slide 19）：

       A   C   G   T
   A   2  -7  -5  -7
   C  -7   2  -7  -5
   G  -5  -7   2  -7
   T  -7  -5  -7   2

2. Gap Penalty（缺口惩罚）
   插入（insertion）或缺失（deletion）称为 indel，每出现一次 gap，要扣分。

   常用模型：Affine Gap Penalty

   $$
   Penalty = d + (n-1)×ε
   $$

   • $d$：gap opening penalty（开缺口惩罚）
   • $ε$：gap extension penalty（延伸惩罚）

   这样可以鼓励连续的缺口，而不是分散多个小缺口。

最终总得分为：

$$
Final\ Score = Σ(substitution\ scores) - Σ(gap\ penalties)
$$

---

Part 2. Global & Local Alignment — Dynamic Programming

1. Pairwise Sequence Alignment: Mathematical Formulation

Input: Two sequences S₁ and S₂
Parameter: Scoring function f for substitution & gaps
Output: The alignment ali(S₁, S₂) that gives maximal score

$$
\max_{ali} f(ali(S_1, S_2))
$$

但暴力穷举是不可能的。
因为两条长度为 n 的序列，其所有可能排列方式为：

$$
\binom{2n}{n} = \frac{(2n)!}{(n!)^2}
$$

当 n = 300 时，这个数约为 $7×10^{88}$，远超可计算范围。
因此我们需要一种高效算法——Dynamic Programming（动态规划）。

---

2. DP Recurrence Formula (for Global Alignment)

假设我们要对齐两个序列：

• x = x₁x₂…xₘ
• y = y₁y₂…yₙ

定义：

• $F(i,j)$：前 i 个字符的 x 与前 j 个字符的 y 的最佳得分。
• $s(x_i,y_j)$：匹配得分函数。
• $d$：线性 gap penalty。

公式为：

$$
F(i,j) = \max
\begin{cases}
F(i-1,j-1) + s(x_i, y_j) & (x_i \leftrightarrow y_j) \
F(i-1,j) + d & (x_i \leftrightarrow gap) \
F(i,j-1) + d & (y_j \leftrightarrow gap)
\end{cases}
$$

初始条件：$F(0,0)=0$

这个公式在 slide 17–18 展示了详细图示：
箭头表示三种可能的转移方向。

---

3. Example: Global Alignment of AAG vs AGC

（slide 20–25）

1. Substitution matrix:

       A   C   G   T
   A   2  -7  -5  -7
   C  -7   2  -7  -5
   G  -5  -7   2  -7
   T  -7  -5  -7   2

2. Gap penalty d = -5

逐步计算 $F(i,j)$，得到最终最优得分为 -6，
并通过 traceback 得到比对结果：

A A G -
- A G C

这就是一个典型的 global alignment。

---

4. From Global to Local Alignment

Global alignment 适合长度相似的序列，但当我们只关心一部分（如基因中的保守片段）时，需要 local alignment。

Smith–Waterman algorithm (1981)
修改了上面的 DP 公式：
多加一个选项 0，防止负分累积。

$$
F(i,j) = \max
\begin{cases}
F(i-1,j-1)+s(x_i,y_j) \
F(i-1,j)+d \
F(i,j-1)+d \
0
\end{cases}
$$

这样，当某个比对片段的分数掉到 0 以下时，就“重新开始”计算。

Traceback 从矩阵中得分最高的格子开始，
往回走到第一个为 0 的格子，得到局部相似区。

---

5. Global vs Local Comparison

Feature	Global Alignment	Local Alignment
Algorithm	Needleman–Wunsch	Smith–Waterman
Start of Traceback	Bottom-right corner	Highest score cell
End condition	到达矩阵左上角	分数降为 0
Scope	Whole sequences	Local region only
Application	Homologous sequences	Conserved motifs / domains

---

Part 3. NGS Reads Mapping 高通量测序读段比对

1. What is NGS Reads Mapping?

Definition:
NGS (Next-Generation Sequencing) 生成海量短序列（short reads），
reads mapping 是指将这些短序列准确地定位（map）回参考基因组（reference genome）的过程。

典型的输入与输出如下：

项目	内容
Input	Reference genome (数亿碱基) + Millions of reads (短序列, 36–150 bp)
Output	每个 read 在参考基因组中的最佳匹配位置

---

2. Key Challenges 难点

1. 数据量极大 (Data volume)
   人类基因组 ~3×10⁹ bp，而 NGS 输出可达数亿 reads。
   → 暴力比对不可行。

2. 错误与突变 (Errors & Variations)
   每个 read 可能含有测序错误、SNP、indel。
   → 必须允许部分错配（mismatch allowance）。

3. 重复区域 (Repetitive regions)
   染色体中有大量重复序列，
   → 同一个 read 可能匹配多个位置（multi-mapping）。

因此，NGS mapping 的核心任务是：

找到每个 read 在 reference 上最可能的匹配区域，同时保证速度和准确度。

---

3. From Alignment to Mapping

在概念上，mapping 是 alignment 的特例。
但区别在于：

比较项目	Alignment	Mapping
输入	2 条序列	一个 read 与整个参考基因组
目标	计算相似性	确定 read 的位置
复杂度	O(n²)	必须更快（近似 O(n)）
典型算法	Needleman–Wunsch / Smith–Waterman	BLAST / BWA / Bowtie

因此 mapping 算法的策略是：

用“快速近似 + 精确局部修正”的方式代替完全动态规划。

---

4. Strategy Overview 策略概览

在 slide 45–47 中提到 NGS Mapping 的三步法：

1. Seed searching 种子搜索
   先找一段短的精确匹配（称为 seed 或 word）。
   种子长度通常 10–15 bp。

2. Extension 扩展
   从种子开始向两边延伸（允许错配）。

3. Evaluation 评分
   计算匹配得分与统计显著性（E-value / MAPQ）。

这种方法就是著名的 BLAST 思想：Seed–Extend–Evaluate。

---

5. Seed Searching Strategies

在 slide 50 提到两种主要策略：

1. Hashing (哈希法)

   • 建立一个索引表（hash table），
     把 reference 中所有 k-mer 存进去。
   • 对 query 的每个 k-mer 查询是否存在相同键。
   • 可以允许少量错配（通过模糊哈希或扩展种子）。

   优点：实现简单；适合较小数据库。
   缺点：内存占用大；对大基因组效率低。

2. BWT (Burrows–Wheeler Transform)

   • 一种基于后缀数组（suffix array）的压缩索引。
   • 可以快速地在压缩数据中搜索字符串。
   • 是现代 NGS mapper（如 BWA、Bowtie）的基础。

---

Part 4. BLAST Algorithm and Its Variants

1. What is BLAST?

BLAST = Basic Local Alignment Search Tool
它是目前最广泛使用的序列比对程序。
基本思想源自 Smith–Waterman 局部比对，但通过“种子 + 扩展”来大幅加速。

BLAST 解决的问题：

给定一个 query 序列，快速找出数据库中与它最相似的序列。

常见变种：

工具	说明
blastn	核酸序列比对
blastp	蛋白质序列比对
blastx	将 query 的核酸序列翻译后与蛋白数据库比对
tblastn	将蛋白序列与翻译后的核酸数据库比对
tblastx	双向翻译比对（最慢）

---

2. BLAST Workflow

（详见 slide 47–50）

1. Word generation（词表生成）
   将 query 序列拆分成固定长度的短词（w-mer）。

   • 核酸：w ≈ 11
   • 蛋白：w ≈ 3

2. Word list expansion（词表扩展）
   对每个词，找出得分较高的相似词（以允许错配）。

3. Database scanning（数据库扫描）
   在数据库中找到完全匹配的词（即 seed）。

4. Ungapped extension（无缺口扩展）
   从匹配种子向两边延伸，直到得分下降超过阈值。

5. Gapped extension（带缺口扩展）
   对得分较高的候选区域进行更精细的比对（允许 gap）。

6. Evaluation and Output

   • 计算每个匹配的得分与统计显著性（E-value）。
   • 输出最高得分的比对结果。

---

3. Statistical Significance — E-value

E-value 表示“在随机序列中出现该得分或更高得分的预期次数”。

公式：

$$
E = K \times m \times n \times e^{-λS}
$$

其中：

• $S$：alignment score
• $m, n$：序列长度
• $K, λ$：经验常数

直观理解：

• E 越小 → 匹配越显著
• 常用阈值：E < 1e-5 被认为是显著相似

---

4. From BLAST to BWA — BWT Indexing

现代 NGS mapping 工具（如 BWA, Bowtie）
并不使用哈希，而是使用 BWT (Burrows–Wheeler Transform) 索引。

BWT 的核心思想（slide 68–70）：

1. 把序列的所有后缀（suffixes）生成并排序。
2. 取每个后缀前一个字符，形成 BWT 字符串。
3. 这种表示法可以高效压缩，并快速进行子串搜索。

示例（slide 69–70: sequence = “AGGAGC”）：

• 构建所有后缀并排序 → 得到 F 列（First）和 L 列（Last）
• 通过 “L → F mapping” 可以逆推出原始序列。

BWA 在此基础上进一步实现：

通过二分查找在压缩索引中定位所有匹配的 reads，再局部比对验证。

---

5. Mapping Quality (MAPQ)

每个 read 通常有多个可能比对位置。
Mapping quality (MAPQ) 用于衡量当前比对的可信度：

$$
MAPQ = -10 \log_{10}(P_{error})
$$

其中 $P_{error}$ 是该 read 被错误比对的概率。

举例：

• 如果正确比对概率 0.99 → MAPQ = 20
• 如果正确比对概率 0.999 → MAPQ = 30

高 MAPQ 表示该位置更可靠。

---

6. Implementations — In R and Python

slide 71–72 给出了两个实现示例：

(1) R Implementation using Biostrings

在 R 中，可以用 Biostrings 包进行简化版 BLAST：

library(Biostrings)
query <- DNAString("ATCGATCGATCG")
db <- DNAString("CGATCGATCGATC")

match <- matchPattern(query, db)
score <- score(match)

代码中展示了如何：

• 分割成 k-mer
• 匹配并计算 similarity score
• 计算 E-value 并输出排序结果

(2) Python Implementation using Biopython

from Bio.Blast import NCBIWWW, NCBIXML
result_handle = NCBIWWW.qblast("blastn", "nt", "ACGTGAGGCTAGC...")
blast_record = NCBIXML.read(result_handle)
for alignment in blast_record.alignments:
    for hsp in alignment.hsps:
        print(alignment.title, hsp.expect)

该脚本演示了如何通过 NCBI 在线数据库运行 BLAST 查询并解析结果。

---

Part 5. Summary 总结

概念	英文术语	核心思想
Sequence Alignment	序列比对	用动态规划计算相似度
Global Alignment	全局比对	对齐整个序列
Local Alignment	局部比对	找最相似片段
Gap Penalty	缺口惩罚	惩罚 indel 事件
DP Algorithm	动态规划	拆分子问题获得最优解
Seed-Extend	种子扩展策略	BLAST 的核心思想
BWT	Burrows–Wheeler Transform	高效索引与压缩算法
MAPQ	Mapping Quality	比对可信度评分
E-value	Expectation Value	随机匹配显著性评估

---