生信工具

查询蛋白信息使用的数据库

NCBI

NCBI–protein–直接查名字，然后可以下载各种格式的文件 fasta

Uniprot

uniprotKB，有review和unreview之分，review是黄色文件，已经有验证的，unreview是灰色没有验证，很多都是预测的 也可以得到fasta文件

Pfam数据库介绍及使用

Pfam: Home pagepfam.xfam.org

Pfam is a large collection of protein families, represented by multiple sequence alignments and hidden Markov models (HMMs)

蛋白质家族数据库

Domain一个基因转录的蛋白质分子中可以包含多个结构特异并且功能不同的区域，这些区域称之为domain，domain 可以看作蛋白质功能的基本单位，蛋白质的功能由包含的多个domain共同决定，研究domain, 可以更好的研究蛋白质功能，而具有相同结构域的基因往往形成一个基因家族。

基因家族是来源于同一个祖先，由一个基因通过基因重复而产生两个或更多的拷贝而构成的一组基因，它们在结构和功能上具有明显的相似性，编码相似的蛋白质产物。

做Alignment

NCBI–Protein–blastp

但是blast在线可以做的序列长度有限，所以可以选择脱机本地跑，也可以选择其他方法

Clustal Omega

Clustal Omega：在线序列alignment，可以导出多种格式，比如.sto，方便后续hmmer计算

ClustalW和X都是本地的，可以下载（MacOS新版貌似不太行。。）

做完alignment，就知道序列们相似在哪不相似在哪，根据这个特征，HMMER可以将alignment输出的结果建立模型，作为模型工具用于更多序列的检索

这玩意只能匹配两条或更多条序列，也就是说输入的内容必须是2+的序列，只输入一条序列会报错。但是如果我们需要sto文件，该怎么办呢？可以使用seqmagick软件，用来格式转换

pip install seqmagick   #安装
seqmagick convert a.fasta a.sto

蛋白预测

Prodigal

prodigal是一个你给了他核酸序列，他能给你翻译成氨基酸序列，然后把这个序列也进行alignment和hmmer的一顿操作的一个工具

为啥他可以这么搞呢？

Prodigal is a protein-coding gene prediction software tool for bacterial and archaeal genomes. The acronym stands for PROkaryotic DYnamic Programming Genefinding ALgorithm. Dictionary.com provides several definitions of the word “prodigal”.

这是一个针对于原核生物的工具，因为原核生物没有内含子，所以直接翻译就ok，但是真核生物因为内含子的存在，就得按照不同物种进行预测，很麻烦

prodigal -i inputfile -a outputfile

伪代码

• 读入序列
• 定位基因组中所有开始和终止位点
• 扫描所有阅读框，记录开放阅读框中G、C的个数
• 构建基于开放阅读框ORF长度和GC位置的表格偏见模型
• 记录每个模型起始子和终止子overlap小于60bp的最高分
• 做第一次动态编程，连接分数偏见的模型，连接节点
• ……

HMMER: biosequence analysis using profile hidden Markov models

http://hmmer.org/

搜索同源序列，进行序列比对。

HMMER用于在序列数据库中搜索序列同源物，并进行序列比对。它使用被称为轮廓隐马尔可夫模型(profile hmm)的概率模型实现方法，搜索远程同源物。

这里我们用于将预测出的蛋白进行hmmer运算，从而预测得到更多的蛋白。

hmmbuild

build a profile HMM from an alignment

Synopsis

hmmbuild [options] hmmfile alignfile

Description

hmmbuild reads a multiple sequence alignment file alignfile , builds a new profile HMM, and saves the HMM in hmmfile.

alignfile may be in ClustalW, GCG MSF, or SELEX alignment format.

By default, the model is configured to find one or more nonoverlapping alignments to the complete model. This is analogous to the behavior of the hmmls program of HMMER 1. To configure the model for a single global alignment, use the -g option; to configure the model for multiple local alignments a la the old program hmmfs, use the -f option; and to configure the model for a single local alignment (a la standard Smith/Waterman, or the old hmmsw program), use the -s option.

hmmsearch

hmmsearch hmmfile.hmm faalist > outputfile

序列处理

seqtk subseq

序列提取，按照namelist中的序列从fastqfile.fa中进行提取

seqtk subseq fastqfile.fa namelist > output.fa 

BWA

bwa mem

# 根据reference genome data(e.g. ref.fa) 建立 Index File：
bwa index ref.fa
# 输出的sam文件还包含了双端比对的相关信息
bwa mem $REF $R1 $R2 > bwa.sam
# 在这里，$1和$2可以是fastq格式，也可以是压缩文件，可以说是很人性化了，但是生成的文件非常大，一定要注意好存储空间

Bowtie2

samtools

先把mapping上的序列提取出来，然后把提取的sam文件压缩成bam文件，以节省空间

IGV

查看mapping结果