Bioinformatics Data

约 2860 个字 67 行代码预计阅读时间 11 分钟

在生信分析当中，我们要面对非常大的数据量（对于最新的各种测序方法来说，百万甚至千万的数据点都不算什么），这种数据量大到我们不能以传统的方式去处理这些数据，比如人手去整理（诗人握持）。

概括来说，处理生物信息学的数据有以下几个难点：

检索数据：如何获取我们需要的数据？有时候我们甚至要获取很多数据。
数据完整性：数据可能在传输的任何阶段出现损坏，我们如何保证拿到的数据是原来的样子？
数据压缩：生信所面对的数据量太大了，为了性能考虑，我们该如何压缩数据？

Retrieving Bioinformatics Data

本节将介绍在Unix shell里常用的两个命令行工具：wget和curl。

Downloading Data with wget and curl

wget

wget是一个强大的命令行下载工具，支持HTTP/HTTPS、FTP等等协议，很多发行版都自带了它，当然也可以通过包管理器安装：

❯ sudo apt-get update
❯ sudo apt-get install wget

用wget下载基因组文件：

❯ wget http://hgdownload.soe.ucsc.edu/goldenPath/hg19/chromosomes/chr22.fa.gz
--2025-05-27 21:29:33--  http://hgdownload.soe.ucsc.edu/goldenPath/hg19/chromosomes/chr22.fa.gz
Resolving hgdownload.soe.ucsc.edu (hgdownload.soe.ucsc.edu)... 128.114.119.163
Connecting to hgdownload.soe.ucsc.edu (hgdownload.soe.ucsc.edu)|128.114.119.163|:80... connected.
HTTP request sent, awaiting response... 200 OK
Length: 11327826 (11M) [application/x-gzip]
Saving to: ‘chr22.fa.gz.1’

chr22.fa.gz.1                 100%[=================================================>]  10.80M  3.53MB/s    in 3.1s

2025-05-27 21:29:37 (3.53 MB/s) - ‘chr22.fa.gz.1’ saved [11327826/11327826]

如何过程顺利，你可以在当前工作目录下发现文件chr22.fa.gz，使用你喜欢的解压工具解压它（我使用zsh，所以用的x），就可以获得包含基因组信息的文件了。（应该是chr22.fa）

在下载大文件的时候，FTP的表现会比HTTP更加优秀，这也是UCSC所推荐的方式。

wget还可以完成简单的FTP或HTTP验证，只需要使用其的--user和--ask-password即可。

wget的另一个强大之处在于，其支持递归的下载页面文件，就像一个爬虫一样，这需要我们启用--recursive参数。而为了防止wget无限的下载内容，我们可以通过--level参数项控制递归深度，默认为5。

下面我们从实例来看看wget的递归模式：

❯ wget --accept "*.aso.gz" --no-directories --recursive --no-parent \
https://download2.cncb.ac.cn/genbase/GenBank/daily/2024/

使用以上命令可以直接下载国家基因组科学数据中心提供的所有GenBank基因组信息，但是注意，是所有！！！

以上命令中，我们指定了--no-directories和--no-parent来阻止wget下载文件夹和父目录内容。同时使用--accept限定了目标文件。

当然，当wget递归执行时候，其会在短时间内对远程服务器做出大量请求，给服务器造成压力，并有可能导致我们的IP被封禁。为了避免这种情况，我们可以通过--limit-rate选项来控制wget的下载速度

curl

curl与wget行为类似，但是不同在其默认将文件写入标准输出，所以要像wget一样下载文件，我们可以通过curl的--output参数指定文件名，或者将输出重定向至某个文件内：

❯ curl http://hgdownload.soe.ucsc.edu/goldenPath/hg19/chromosomes/chr22.fa.gz > chr22.fa.gz
  % Total    % Received % Xferd  Average Speed   Time    Time     Time  Current
                                 Dload  Upload   Total   Spent    Left  Speed
100 10.8M  100 10.8M    0     0  3191k      0  0:00:03  0:00:03 --:--:-- 3191k

比起wget，curl的优势在其支持更多的协议，比如SFTP和SCP，如果启用了-L/--location，curl可以自动处理页面重定向，直到指定的下载链接。并且，curl本身是一个库，它的功能已经被RCurl和pycurl所封装。

Rsync and SCP

wget和curl比较适合从网络上下载对应的文件，而当我们需要在多个设备之间同步某些文件的时候，rsync是一个不错的选择。

rsync在某种意义上很像Git，它确实可以管理文件的版本，并且可以操作文件，是一种理想的归档工具。

rsync基本的语法格式为：

❯ rsync [OPTION]... SRC [SRC]... DEST

假如我们有一个远程主机，其的IP地址为172.26.50.20，我们想将本地的一个文件同步到远程主机的用户yangshu233上去，命令如下：

❯ rsync -avz -e ssh /path/to/my/file/ yangshu233@172.26.50.20:/home/abc/data

特别地，对于rsync来说，源目录SRC的路径最后是否存在斜杠/的行为是区别的，当存在斜杆的时候，比如project/data/，rsync只会将data/下的所有内容复制到指定的远程目录；若不存在斜杠，如project/data，rsync会拷贝整个目录，包括data在内，至指定远程目录下。

当然，有时候我们只是想要传递一些文件，而无需在一段时间内对其进行同步或归档操作，这时候使用rsync有点大材小用，而scp是一个更方便的选择，它通过SSH连接到远程主机，并以类似cp的方式拷贝文件。

scp /my/data/file user@host:/home/user/data

特别地，scp总是拷贝整个源目录，与源目录路径最后是否有斜杠/无关。

Data Integrity

在平常生活中，我们可能很少感受到在网络数据传输中出现的损坏问题，很多时候是因为我们所使用的软件可以自动处理这些错误，或者只是单纯的因为传输的数据量小而短，出错的概率不大。

但是在生物信息学中，一段基因组的信息，或者其他数据可能非常大，这需要我们花费较多时间下载它们，这大大增加了因为一些莫名其妙的因素而导致数据损坏的可能。为了分析工作的顺利开展，我们有必要使用一些手段来确定我们拿到的数据是正确的，这需要校验和这个手段。

Note

我们通过某一种算法，对文件的内容进行压缩和编码，所得出的有关于该文件的唯一标识符，就是该文件的校验和。对于每个文件来说，校验和都是唯一的。即使只有一个比特的改变，校验和也会完全不同。

校验和不仅可以用于对数据完整性的保证，还可以参与到数据的版本控制当中。比如我们的程序需要针对不同版本的数据生成不同的结果，我们该如何确定两组数据是否是同一个版本？该如何将结果和数据联系在一起？这些都可以通过校验解决！

SHA and MD5 Checksums

SHA和MD5是两种常见的校验和算法，其中SHA-1被用于Git生成提交ID，而MD5则是一个较老但广受欢迎的算法，如果一个网站提供了某个文件的校验和，那么它很可能是由MD5生成的。

让我们来看看SHA-1校验和算法：

❯ echo "hello world" | shasum -a 1
22596363b3de40b06f981fb85d82312e8c0ed511  -
❯ echo "hello worl" | shasum -a 1
10bbef0200dced22b2c707c8621a81d571d224f0  -
❯ echo "hello world" | shasum -a 1
22596363b3de40b06f981fb85d82312e8c0ed511  -

shasum是生成SHA校验和的命令，这里我用-a参数指定了使用SHA-1，我们可以发现，前后两个字符串尽管只差了结尾的一个d，生成的校验和也截然不同，同时无论何时，只要使用相同的算法，SHA-1对相同的字符串生成的校验和是完全相同的。这意味着其具有稳定性和唯一性，是一个非常可靠的校验和算法。

我们还可以将文件传给shasum来生成校验和：

❯ shasum -a 1 chr22.fa
4dd1c94043fd40d31e6ddc16f34075b78044986d  chr22.fa

当然，如果我们有大量的文件需要进行进行校验，手动处理会非常麻烦。对此，shasum提供了通过checksum file进行校验的方法：

# 生成checksum file
❯ shasum data/*fastq > fastq_checksums.sha
# 使用checksum file进行校验
❯ shasum fastq_checksums.sha

在Unix shell中，同样提供了md5sum来计算MD5的校验和，在一些老旧的服务器上，可能使用的是chsum或者sum。

Looking Differences Between Data

校验和只告诉我们两个文件是否不同，但是不能告诉我们文件内容的具体差异是什么，在Unix shell中，我们有diff这命令，它可以高效地比较两个文件之间的差异，并以“块”为单位输出所有的不同。

❯ diff -u gene-1.bed gene-2.bed
--- gene-1.bed  2025-05-28 16:20:59.365511235 +0800
+++ gene-2.bed  2025-05-28 16:20:59.365511235 +0800
@@ -1,22 +1,19 @@
 1      6206197 6206270 GENE00000025907
 1      6223599 6223745 GENE00000025907
 1      6227940 6228049 GENE00000025907
+1      6222341 6228319 GENE00000025907
 1      6229959 6230073 GENE00000025907
-1      6230003 6230005 GENE00000025907
 1      6233961 6234087 GENE00000025907
 1      6234229 6234311 GENE00000025907
 1      6206227 6206270 GENE00000025907
 1      6227940 6228049 GENE00000025907
 1      6229959 6230073 GENE00000025907
-1      6230003 6230073 GENE00000025907
+1      6230133 6230191 GENE00000025907
 1      6233961 6234087 GENE00000025907
 1      6234229 6234399 GENE00000025907
 1      6238262 6238384 GENE00000025907
-1      6214645 6214957 GENE00000025907
 1      6227940 6228049 GENE00000025907
 1      6229959 6230073 GENE00000025907
-1      6230003 6230073 GENE00000025907
 1      6233961 6234087 GENE00000025907
 1      6234229 6234399 GENE00000025907
-1      6238262 6238464 GENE00000025907
 1      6239952 6240378 GENE00000025907

此处我使用-u，使得输出的格式类似于git diff，它可以提供更多的上下文信息。我们简短了解一下其输出内容的含义：

头两行给出了两个文件的基本信息，其中---代表原文件，+++代表修改后的文件。
@@...@@之间的内容代表原文件和修改后文件的修改起始位置。
接下来diff给出了具体的差异之处，其中1代表相同，-1代表删除，+1代表新增。如果-1行后紧跟+1行，通常意味着这行被修改。

特别地，diff的输出可以被重定向到文件当中，这个文件被我们称作补丁，其作为修改一个纯文本的指令，在Unix shell提供了patch命令，可以根据补丁来应用更改。

Compressing Data and Working with it

某些大型生物信息学的数据集非常之大，可以达到数十GB，如果它们未经任何处理，将会占用大量的磁盘空间，这可能会严重影响磁盘性能，或者干脆使得某些需要一些缓存空间的程序没法运行。

于是，我们有时候必须对本地的数据进行压缩，幸运的是，不仅Unix shell提供了gzip这样的压缩工具，cat，grep和less也有支持压缩数据的变体，而Python更是提供了gzip模块，支持对压缩数据的读写。

gzip

Unix shell中常用的两种压缩工具是gzip和bzip2，前者压缩和解压速度更快，后者的压缩率更高。鉴于以上特点，我们一般选择gzip来压缩我们的数据，而bzip2更多用于长期数据的压缩存储。

gzip可以从标准输入获取数据并压缩，比如我们有一个程序输出序列数据，我们可以这么做：

❯ program data.fastq | gzip > out_data.fastq.gz

gzip还可以直接压缩文件：

❯ gzip data.fastq

此时你的工作目录下应该会出现一个data.fastaq.gz文件，再次使用gzip可以将其解压：

❯ gunzip data.fastaq.gz

这样你可以获得文件data.fastaq。

Note

gzip和gunzip在压缩和解压文件的时候，会将原压缩/解压缩文件删除。

当然，可以通过gzip的-c参数将压缩后的文件输出到标准输出中，由我们手动重定向至特定文件中，这不仅可以保留原文件，还可以实现向已压缩的文件追加压缩内容的功能。

❯ gzip -c data.fastaq > data.fastaq.gz
❯ gzip -c data2.fastaq >> data.fastaq.gz

Working with Compressed Data

在Unix shell中，grep的变体zgrep，cat的变体zcat和less变体zless都可以直接处理压缩文件，并且其用法与原命令非常类似，参数也是一致的。

同时，很多生物信息学工具都支持压缩文件的操作，如果某个程序或工具不支持直接读取压缩数据，可以通过zcat配合管道符传递输入。

由于在操作数据的时候需要先解压缩，这可能会产生一些CPU性能上的损失，但是总体上利大于弊，因为大量的磁盘空间被节省，而CPU的消耗几乎可以忽略不计。

Reproducibly Downloading Data

很多时候，我们的工作不仅局限于当下，今天我们下载的数据，在一段时间后可能需要再次下载来分析，或者和我们合作的其他人也需要下载这个数据。我想，如果能详细记录下载的必要信息和文件的校验和，将对将来的研究很有帮助。

这份Markdown至少应该包含有：

下载文件名，下载源和时间，同时附上下载链接或者下载命令。
相应文件的校验和。