已知PMID，得到参考文献的格式
最近太忙了，抽空写了个程序，可以批量从pmid得到参考文献。
pmidlist = open('pmid.txt', 'r')
output = open('results.txt', 'w')
def pmid2ref(pmid):
import requests
from bs4 import BeautifulSoup
html = requests.get('https://www.ncbi.nlm.nih.gov/pubmed/' + str(pmid) + '/')
soup = BeautifulSoup(html.text, 'lxml')
title = soup.title.string.split('- PubMed')[0]
info = so...
GATK推荐的germline call snp+indel流程
GATK4推荐流程。
其实还没完全搞懂。
不过把变异call出来是没问题的！
没有和其他软件对比过正确率。。。
#1 首先把原始数据处理成可以用 的bam
参考推荐的数据准备流程。
#2 同样的，配置好软件和环境
像这样。
sample=sample.pre.analysis.bam
gatk=/software/gatk-4.0.0.0/gatk
reference=/database/GATK/hg19/ucsc.hg19.fasta
omni=/database/GATK/hg19/1000G_omni2.5.hg19.vcf
indel=/database/GATK/hg19/1000G_phase1.indels.hg19.vcf
dbsnp=/database/G...
GATK推荐的数据预处理流程
GATK4的推荐流程哦。这个流程是从ubam文件或fastq文件得到可以用于找变异的bam文件。
以后GATK大概会把bwa和picard完全内嵌。
下面的命令其实还有很多参数的，不过一般来说用这些参数就够。
原始数据要用经过了质控之后的数据！
放两个用作测试学习的原始数据。
点这里下载！
#1 先把需要用到的软件还有数据库下载好
首先是软件。当然是要用最新的啦，跟上时代。
GATK4.0
picard2.17.4
bwa0.7.17
然后是数据库。
用的是GATK提供的hg19。
地址在这：点击进入
#2 把软件都加到环境中，或者在shell脚本中设定好
像这样。
reference=/database/hg19/ucsc.hg19.fasta
indel1=/d...
从TCGA中获得甲状腺癌的相关突变
很忙很累还更新的我。
TCGA是一个癌症相关的数据库。可以弄到很多癌症的基因数据样本。
点点看TCGA。
进去之后，点右边那个蓝色的Launch Data Portal
可以看到样本量还是挺多的。
我们可以点右边的那个人的甲状腺（Thyroid）。
就进入了一个统计界面。
可以看到截至目前的相关突变位点有11128个。
单击Mutations，再点击JSON就可以把所有的突变位点下载下来。
然后我用下面这个脚本，提取出了所有的位点。
inputfile = open('Thyroid_mutation.json', 'r')
outputfile = open('results.txt', 'w')
for line in inputfile:
if line.sta...
多序列比对软件mafft
2018年的第一次更新。
mafft是一个用来做比对的软件。最近我在研究HLA的分型（用的一代）的时候想用来试试。
先看看官网。
可以看到目前的最新版本是7.313。
我装的是这个，linux版本的。
官网提供了按照的方法：
# 在root下
rpm -Uvh mafft-xxxxx-rpm
exit
# if necessary
rehash
装好之后可以用:
mafft -help
查看帮助文档。
最简单的使用方式就是：
mafft input > output
比如说我有一个HLA-C区的fasta文件，点这里下载。
然后运行下面这个就可以了。
mafft C_nuc.fasta > output.txt
建立阿尔兹海默病的panel
简单的，做panel的，例子
panel就是关于某个疾病的位点，一般包括这个疾病的致病基因位点，药敏基因位点，毒副作用基因位点，风险基因位点这些。
很多时候都要靠手工去查询文献来手动添加。或者，可以去购买商业化的疾病panel。
一个panel，最重要的是可信。可信的意思是，位点要有相关文献的支持。
然后位点要的是参考基因组上的绝对位置。
例如hg19的chr1:0000001这样的。
下面利用这个公开的老年痴呆数据库来建立一个panel。
数据库
这个数据库的有点在于，每个位点都有pmid，而且标注了位点是致病的还是良性的。
首先，我们通过浏览器右键查看源代码简单粗暴的把这个网页弄下来。
当然，页面很多的情况下是写一个爬虫程序来把网页弄下来，但是我们需要的只是一个页面，...
使用Transvar来找位点
很多时候，看文献啊报告这类东西的时候作者因为某些原因，不会给出明确位点，只会给基因名和碱基突变或蛋白质突变。例如，给出 ITPA c.94C>A 这样的突变。
这里介绍一个可以凭这点信息找到绝对位置（染色体+物理长度）的工具：Transvar
Transvar是一个用python写的开源软件，我们可以下载下来使用，但是下载的话还要下载数据库比较麻烦。所以可以选择使用网页版。
Transvar Web。
网页版还是一个做了用户界面的，非常贴心而且非常容易使用。看一眼就会的那种。
如上图，Select a task可以选择protein（aachange）、cDNA、gDNA等。
Select a reference genome可以选择参考基因组。
Select on...
建立阿尔兹海默病的数据库
要弄一个项目，首先要找到相关的panel。比如说，要做老年痴呆的项目，首先就应该去找关于阿尔兹海默病和额颞叶痴呆的相关基因位点。
刚好的是，发现一个统计了相关文献研究的网站。
这里molgen.vib-ua.be。
然后选择基因可以看到这样的一个界面。
首选的弄下来的方法，当然是复制粘贴。
然而，复制粘贴格式会乱掉。
只好用脚本处理了。
我们可以找到这个网页的源代码。
点这里
下面是我用来处理的python脚本：
from bs4 import BeautifulSoup
inputFile = open('AD.html', 'r')
outputFile = open('result.txt', 'w')
soup = BeautifulSoup(input...
共计 232 篇文章，29 页。