基因组测序分析linux,基因组数据的重测序分析
我通过查资料获得已知达松维尔拟诺卡氏菌亚种(cardiopsis dassonvillei subsp. dassonvillei)的基因组原始测序序列,我想知道这个亚种与达松维尔拟诺卡氏菌(Nocardiopsis dassonvillei )的基因组相比有哪些不同。
1.需要的软件
软件名:Aspera 版本号:3.6.2.117442
软件名:sratoolkit 版本号:2.9.2
软件名:FastQC 版本号:0.11.7
软件名:Trimmomatic版本号:0.38
软件名:bwa 版本号:0.7.17-r1188
软件名:samtools 版本号:1.7
软件名:Annovar 版本:$Date: 2017-07-17 01:16:48 -0400 (Mon, 17 Jul 2017) $
2.数据下载
~/.aspera/connect/bin/ascp -T -i /home/lizeguo/.aspera/connect/etc/asperaweb_id_dsa.openssh -k 1 -l 200m anonftp@ftp-private.ncbi.nlm.nih.gov:/sra/sra-instant/reads/ByRun/sra/SRR/SRR022/SRR022534/SRR022534.sra ./Seqs/
~/.aspera/connect/bin/ascp -T -i /home/lizeguo/.aspera/connect/etc/asperaweb_id_dsa.openssh -k 1 -l 200m anonftp@ftp-private.ncbi.nlm.nih.gov:/genomes/all/GCA/001/877/055/GCA_001877055.1_ASM187705v1/GCA_001877055.1_ASM187705v1_genomic.fna.gz ~/Seqs/
基因组gff文件
cd ~/Seqs
wget ftp://ftp.ncbi.nlm.nih.gov/genomes/all/GCA/001/877/055/GCA_001877055.1_ASM187705v1/GCA_001877055.1_ASM187705v1_genomic.gff.gz
3.主要分析步骤和结果
(一)序列与参考基因组的比对
1.数据文件的格式转换
在NCBI的SRA数据库检索序列号SRR022534查看测序数据类型为454的单末端测序
2018-12-20 10-57-15屏幕截图.png
fastq-dump SRR022534.sra
2018-12-20 11-04-11屏幕截图.png
2.质量评估
mkdir fastqc_result
fastqc SRR022534.fastq
mv *.zip *.html fastqc_result/
2018-12-20 11-16-21屏幕截图.png
3.测序数据的数据过滤
mkdir trim_out
java -jar ~/BioSofts/Trimmomatic-0.38/trimmomatic-0.38.jar SE -phred33 SRR022534.fastq ./trim_out/SRR011534_out.fastq.gz ILLUMINACLIP:/home/lizeguo/BioSofts/Trimmomatic-0.38/adapters/TruSeq2-PE.fa:2:30:10 SLIDINGWINDOW:5:20 LEADING:20 TRAILING:20 MINLEN:75
2018-12-20 11-29-33屏幕截图.png
4.建立参考基因组索引
gunzip GCA_001877055.1_ASM187705v1_genomic.fna.gz
bwa index GCA_001877055.1_ASM187705v1_genomic.fna
2018-12-20 11-37-59屏幕截图.png
5.测序数据比对到参考基因组得到sam文件
bwa mem GCA_001877055.1_ASM187705v1_genomic.fna trim_out/SRR011534_out.fastq.gz >bwa_mem_SRR011534.sam
6.sam文件转换为bam文件
samtools faidx GCA_001877055.1_ASM187705v1_genomic.fna
samtools view -bhS -t GCA_001877055.1_ASM187705v1_genomic.fna.fai -o bwa_mem_SRR011534.bam bwa_mem_SRR011534.sam
7.为bam文件排序
samtools sort bwa_mem_SRR011534.bam -o bwa_mem_SRR011534.sorted.bam
8.为bam文件建立索引
samtools index bwa_mem_SRR011534.sorted.bam
9.显示基因组比对情况
samtools tview bwa_mem_SRR011534.sorted.bam GCA_001877055.1_ASM187705v1_genomic.fna
reads比对情况
10.测试参考基因组每个位点或一段区域的测序深度
samtools depth bwa_mem_SRR011534.sorted.bam >>depth.txt
less depth.txt
每个位点或区域的测序深度
11.统计比对结果
samtools flagstat bwa_mem_SRR011534.sorted.bam
比对结果
结果显示有458257个碱基,有441750个碱基匹配上了,占比96.40%.
(二)变异位点的检测
1.去除PCR重复
samtools rmdup bwa_mem_SRR011534.sorted.bam bwa_mem_SRR011534_nopcr.bam
2.生成bcf文件
samtools mpileup -gf GCA_001877055.1_ASM187705v1_genomic.fna bwa_mem_SRR011534_nopcr.bam >bwa_mem_SRR011534.bcf
3.基因变异检测,得到bwa_mem_SRR011534.snps.bcf文件
bcftools call -vm bwa_mem_SRR011534.bcf -o bwa_mem_SRR011534.variants.bcf
bcftools view -v snps,indels bwa_mem_SRR011534.variants.bcf > bwa_mem_SRR011534.snps.vcf
less bwa_mem_SRR011534.snps.vcf
2018-12-31 13-47-46屏幕截图.png
4.变异位点的过滤
bcftools filter -o bwa_mem_SRR011534.snps.filtered.vcf -i 'QUAL>20 &&DP>5' bwa_mem_SRR011534.snps.vcf
(三)变异基因注释
1.生成annovar输入文件
convert2annovar.pl -format vcf4 bwa_mem_SRR011534.snps.vcf > bwa_mem_SRR011534.snps.avinput
2.自定义注释数据库
gunzip GCA_001877055.1_ASM187705v1_genomic.gff.gz
2.1gff3文件转为GenePred文件
wget http://hgdownload.soe.ucsc.edu/admin/exe/linux.x86_64/gff3ToGenePred
chmod 777 gff3ToGenePred
./gff3ToGenePred -useName GCA_001877055.1_ASM187705v1_genomic.gff 7055-genome_refGene.txt
2.2修改GenePred文件
cut -f 12 7055-genome_refGene.txt >column1.txt
cut -f 2-15 7055-genome_refGene.txt >column_else.txt
paste column1.txt column_else.txt >7055-genome_new_refGene.txt
2.3为建立每个基因与编码序列对应文件
retrieve_seq_from_fasta.pl -format refGene -seqfile GCA_001877055.1_ASM187705v1_genomic.fna -outfile 7055-genome_new_refGeneMrna.fa 7055-genome_new_refGene.txt
2.4拷贝数据库文件到annovar安装目录humandb文件夹
cp 7055-genome_new_refGene* ~/BioSofts/annovar/humandb/
3.注释变异基因位点,生成avinput.variant_function和avinput.exonic_variant_function后缀的两个结果文件
annotate_variation.pl --geneanno --dbtype refGene --buildver 7055-genome_new bwa_mem_SRR011534.snps.avinput ~/BioSofts/annovar/humandb/
2018-12-31 15-13-41屏幕截图.png
4.查看结果
less bwa_mem_SRR011534.snps.avinput.variant_function
2018-12-31 15-15-49屏幕截图.png
less bwa_mem_SRR011534.snps.avinput.exonic_variant_function
2018-12-31 15-15-49屏幕截图.png
基因组测序分析linux,基因组数据的重测序分析相关推荐
- mysql门店分析_用MySQL分析某餐饮业数据 并用杜邦分析图标展示
用MySQL分析某餐饮业数据 展示结果如下: 数据分析 1.餐饮业订单单号详情数据 682条订单信息 2.具体点菜明细数据 ...
- 基因组组装程序linux,基因组组装软件SOAPdenovo安装使用
一.下载并安装这个软件 下载地址进下面,但是下载源码安装总是很困难,我直接下载bin文件可执行程序. 解压进入目录 首先make 然后make install即可 安装总是失败,我也不知道怎么回事,懒 ...
- linux日记的监控与分析,linux下apache日志监控与分析——webalizer与awstat
(一) Webalizer 从webalizer官网(http://www.webalizer.org/)我们可以看到对其做的如下说明,从中我们可以对webalizer有一个简单的了解 The W ...
- linux键盘驱动程序分析,linux设备驱动之键盘驱动分析
一:前言: 在分析intel8042芯片驱动的时候,对其中断处理的后续流程还没有讨论完.在本章以键盘通道为索引讲述intel8042的后续处理,这部份内容实际上是独立的键盘驱动部份,多种型号的键盘都有 ...
- crash工具分析linux内核,如何使用crash工具分析Linux内核崩溃转储文件
满意答案 使用 crash 的先决条件 1. kernel 映像文件 vmlinux 在编译的时候必须指定了 -g 参数,即带有调试信息. 2. 需要有一个内存崩溃转储文件(例如 vmcore),或者 ...
- linux软中断分析,linux操作系统下的软中断问题分析_linux教程
前几天把硬中断完成了,这几天在做软中断的事情.现在对linux2.6.12的软中断机制有了一个基本的认识.在此把它的软中断的一些问题做一个记录. 系统中有一个softirq_action结构,它定义在 ...
- linux io 分析,Linux 下的I/O性能分析 iotop
Linux下的IO统计工具如iostat, nmon等大多数是只能统计到per设备的读写情况, 如果你想知道每个进程是如何使用IO的就比较麻烦. 当然如果你会systemtap, 或者blktrace ...
- linux pinctl 分析,Linux驱动架构之pinctrl子系统分析(一)
1.前言 在嵌入式系统中,许多SoC的内部都包含了pin控制器,通过芯片内部的pin控制器,我们可以配置一个或者一组引脚的状态和功能特性,Linux内核为了统一各SoC厂商的引脚管理,提供了pinct ...
- 全基因组重测序基础及高级分析知识汇总
全基因组重测序基础及高级分析知识汇总 oddxix 已关注 2018.09.20 17:04 字数 11355 阅读 212评论 0喜欢 6 转自:http://www.360doc.com/cont ...
- Crawl:利用bs4和requests爬取了国内顶级某房源平台(2020年7月2日上海二手房)将近30*100多条数据并进行房价分析以及预测
Crawl:利用bs4和requests爬取了国内顶级某房源平台(2020年7月2日上海二手房)将近30*100多条数据并进行房价分析以及预测 目录 利用bs4和requests爬取了国内顶级某房源平 ...
最新文章
- 调查:“AI红娘”来了,有一半人选择相信TA
- 干货!OKR最基础的操作九个原则
- java this context,java – Spring XML中applicationcontext的“this”引用
- HTML5----CSS显示半个字符
- 浅析企业网站建设的基本费用组成
- webui框架的利与弊
- 什么是mysql的游标_数据库中的游标到底是什么意思
- 【.NET进程通信】初探.NET中进程间通信的简单的实现
- github和pypi_如何将GitHub用作PyPi服务器
- win7怎么关uac_win7电脑如何关闭uac用户账户控制
- freeswitch java_freeSWITCH
- [白话解析] Flink的Watermark机制
- smtp协议与pop3协议
- 伊利诺伊香槟分校计算机排名,伊利诺伊大学香槟分校计算机工程排行业界内最整体分析...
- 小白鼠测试---VR头戴设备-暴风魔镜4
- 【ERROR】Error: Error trying to upgrade business network. Error: No valid responses from any peers.
- 森林防火三维地理信息系统
- javax.mail.MessagingException: Could not connect to SMTP host: smtp.163.com, port: 465;
- 计算机动画可分为二维和三维动画,二维动画与三维动画设计的区分
- 数论-FTT 和 NTT