A Byte of Biology

（更新时间：2023-11-26）

物种、分类群

我们的论文在动植物中对HapHiC进行了测试和验证，列表如下：

- 阅读完整内容 -

（更新时间：2023-11-24）

写在前面

近期，我们新开发的基因组挂载工具HapHiC在GitHub上公开了代码，相关论文的预印本也发布在了bioRxiv上，并有幸收到了大家热烈的反馈。然而，在这个过程中我发现虽然GitHub上已经有一份使用说明，许多朋友仍然对HapHiC的使用和原理不太了解。这可能是因为GitHub上的说明过于简单，而HapHiC的设计和参数又较为复杂的缘故。同时，越来越多的中国人加入到了基因组学研究的队伍中，并发表了大量高水平的研究成果。因此，我写下中文版的HapHiC系列教程与方法解读，供大家查阅和参考，希望能对HapHiC的用户有所帮助！

- 阅读完整内容 -

在Python中需要初始化字典的时候我一般会调用collections中的defaultdict来实现，这样可以减少代码行数，并且避免更复杂的逻辑判断，比如：

a_list = ['a', 'b', 'c', 'a', 'd']

# 常规方法
count_dict = {}
for a in a_list:
    if a in count_dict:
        count_dict[a] += 1
    else:
        # 如果dict中不存在该key，则需要赋值来初始化
        count_dict[a] = 1

# 使用defaultdict进行初始化只需要4行代码，且不需要if...else进行逻辑判断
import collections
count_dict = collections.defaultdict(int)
for a in a_list:
    count_dict[a] += 1

- 阅读完整内容 -

一个sratoolkit的报错

想到Linux中~和$HOME的区别这个问题，是因为在集群上安装了sratoolkit，在使用fastq-dump时遇到的一个报错：

$ fastq-dump --split-3 A549-DRR016695.1
2019-11-30T16:49:38 fastq-dump.2.9.6 err: path not found while searching dynamic library within file system module - failed to open SRA manager

=============================================================
An error occurred during processing.
A report was generated into the file '/work/xxxxx/ncbi_error_report.xml'.
If the problem persists, you may consider sending the file
to 'sra-tools@ncbi.nlm.nih.gov' for assistance.
=============================================================

同样fasterq-dump也无法正常使用：

$ fasterq-dump --split-3 A549-DRR016695.1
2019-11-30T16:50:19 fasterq-dump.2.9.6 err: cmn_iter.c cmn_get_acc_type( 'A549-DRR016695.1' ).VDBManagerMakeRead() -> RC(rcFS,rcDylib,rcSearching,rcPath,rcNotFound)

- 阅读完整内容 -

众所周知，描述基因组中位置信息的坐标系通常有1-based和0-based之分。举个例子，要表示某条序列的前5个碱基，如下图的ATGCA，用1-based坐标系就应该是1-5，而0-based坐标系应该是0-5（类似与Python中的切片）:

常见的1-based坐标系的有SAM、GFF/GTF、VCF、WIG、Axt、BLAST比对结果；0-based坐标系的文件有BAM、PSL、MAF（这里指的Multiple Alignment Format）、以及BED等。如果坐标区间范围很大，只做一些简单分析或者可视化，这一个碱基的增减可能影响不大。如果区间很小，比如TSS（Transcription Start Site，转录起始位点）；又或者是要根据编码基因coding region的DNA序列转换成对应蛋白的氨基酸序列，那这一个碱基的差别可能就会带来严重错误。

- 阅读完整内容 -

因为最近博士后入职一堆事情，有一段时间没有更新文章了，内心充满罪恶感，这篇文章也在我草稿中躺了很久了。想起这个主题的起因是前几个月给大名鼎鼎Nature子刊审稿（其实是Scientific Reports），发现文章里做各类相关性分析的时候Pearson和Spearman相关系数随缘混用，基本瞎搞。虽然在这个水平杂志上出现这种问题似乎不会让人感到意外，但是意外的是文章内容做得相当扎实，一作也不是新手，而且还在世界顶尖大学里生物信息学方向的课题组。再到后来入职博士后，发现在我理解里统计学水平应该比生物更强的医学专业学生或者老师，对这个也不是很有概念的样子，所以觉得很有必要拿出来讲讲。

首先分别简单介绍一下我们可能最常用到的Pearson和Spearman相关系数是什么东西，不写公式，尽量说得通俗易懂（其实是怕麻烦）。Pearson相关系数很简单，是用来衡量两个数据集的线性相关程度。而Spearman相关系数不关心两个数据集是否线性相关，而是单调相关，Spearman相关系数也称为等级相关或者秩相关（即rank）。下面几个图看一看应该很容易理解：

- 阅读完整内容 -

本文作为总结篇，将在这里对已经写过的几篇文章内容进行概括，可作为整个系列文章的目录，此外还有一些常见问题的解决方法。加上本篇目前一共写了5篇maftools相关文章，基本上按照整个官方文档使用TCGA-LUAD的数据都跑了一遍，并且加上了一些数据整理、重要参数的解释以及自己遇到的问题，应该可以说是目前最完整的中文资料了。其实maftools本身使用起来很简单，在读入数据之后，基本可视化和数据分析通过1~2行代码就能实现，参考官方文档的example足够了。我写得比较细、花的篇幅较多是因为自己也是初学肿瘤的数据分析，通过学习maftools这样功能丰富的分析工具可以更快速入门。

- 阅读完整内容 -

Maftools系列文章：

1. 《maftools使用方法总结以及常见问题》
2. 《肿瘤变异数据分析和可视化工具maftools：安装和文件格式要求》
3. 《肿瘤变异数据分析和可视化工具maftools：突变数据下载和可视化》
4. 《肿瘤变异数据分析和可视化工具maftools：突变的数据分析》
5. 《肿瘤变异数据分析和可视化工具maftools：CNV的可视化》

上篇文章《肿瘤变异数据分析和可视化工具maftools：突变的数据分析》介绍了使用maftools分析MAF格式的突变数据，但maftools本身功能不限于此，它也可以针对拷贝数变异进行一些分析。

CNV数据下载和处理

还是继续之前TCGA-LUAD的例子。Maftools接受的CNV数据需要是GISTIC输出的结果，包括4个文件——all_lesions.conf_XX.txt、 amp_genes.conf_XX.txt、del_genes.conf_XX.txt以及scores.gistic（其中XX代表置信水平）。但是从TCGA官网能下载到的GISTIC的输出结果只有LUAD.focal_score_by_genes.txt，因此需要下载更上游的DNAcopy的输出结果，然后自己再根据TCGA的流程和参数跑一遍GISTIC。

- 阅读完整内容 -