miRNA生物信息数据分析流程初探(一)
snakemake-qiaseq-mirna流程
在github的snakemake-qiaseq-mirna[1]流程中,MiRNA的分析步骤如下:
- Step 1: Read trimming and UMI detection
- cutadapt[2]
- umi_tools[3]
- Step 2: Read alignment
- Step 3: Deduplication
- Step 4: Read counting
第一步: clone Github的仓库
git clone https://github.com/jounikuj/snakemake-qiaseq-mirna.git
技巧提示
要测试snakemake流程,可以先使用其提供的测试数据。在该仓库中的其测试数据在.test
目录,可以使用snakemake -np --directory .test
命令直接运行,测试流程是否正确,并且查看流程图。
.
├── config.yaml
├── fastq
├── LICENSE.md
├── README.md
├── samples.tsv
├── Snakefile
├── .test
│ ├── config.yaml
│ ├── fastq
│ │ ├── sampleA_S1_L001_R1_001.fastq.gz
│ │ ├── sampleA_S1_L002_R1_001.fastq.gz
│ │ ├── sampleA_S1_L003_R1_001.fastq.gz
│ │ ├── sampleA_S1_L004_R1_001.fastq.gz
│ │ ├── sampleB_S2_L001_R1_001.fastq.gz
│ │ ├── sampleB_S2_L002_R1_001.fastq.gz
│ │ ├── sampleB_S2_L003_R1_001.fastq.gz
│ │ ├── sampleB_S2_L004_R1_001.fastq.gz
│ │ ├── sampleC_S3_L001_R1_001.fastq.gz
│ │ ├── sampleC_S3_L002_R1_001.fastq.gz
│ │ ├── sampleC_S3_L003_R1_001.fastq.gz
│ │ └── sampleC_S3_L004_R1_001.fastq.gz
│ ├── samples.tsv
│ └── .snakemake
│ ├── auxiliary
│ ├── conda
│ ├── conda-archive
│ ├── incomplete
│ ├── locks
│ ├── log
│ ├── metadata
│ ├── shadow
│ └── singularity
└── workflow
├── envs
│ └── env.yaml
├── rules
│ ├── alignment.smk
│ ├── common.smk
│ ├── counting.smk
│ ├── deduplication.smk
│ ├── resources.smk
│ ├── targets.smk
│ └── trimming.smk
├── schemas
│ ├── config.schema.yaml
│ └── samples.schema.yaml
└── scripts
└── merge_read_counts.py
第二步: 查看流程图
在生成流程图的时候,可能需要安装以下软件
sudo apt install graphviz
生成流程图的具体命令是:
snakemake --dag --directory .test | dot -Tpng > dag.png
第三步:用一个样本解释miRNA的流程
使用下面命令生成sampleA
的流程图
snakemake --dag --directory .test results/sampleA/sampleA_counts_dedup.tsv | dot -Tpng > dag_single.png
3.1 数据库下载及建立索引
get_pirna_sequences : piRNA数据下载
curl -s http://www.regulatoryrna.org/database/piRNA/download/archive/v2.0/fasta/hsa.fa.gz | zcat - > resources/pirna.fasta 2> resources/logs/get_pirna_index.log
get_pirna_index:使用bowtie2建立piRNA数据的index
bowtie2-build --threads 1 --quiet resources/pirna.fasta resources/pirna > resources/logs/get_pirna_index.log
get_mirna_sequences : miRNA数据下载
curl -s ftp://mirbase.org/pub/mirbase/CURRENT/mature.fa.gz | zcat - | seqkit grep -r -p ^hsa | seqkit seq --rna2dna > resources/mirna_mature.fasta 2> resources/logs/get_mirna_sequences.log
curl -s ftp://mirbase.org/pub/mirbase/CURRENT/hairpin.fa.gz | zcat - | seqkit grep -r -p ^hsa | seqkit seq --rna2dna > resources/mirna_hairpin.fasta 2> resources/logs/get_mirna_sequences.log
cat resources/mirna_mature.fasta resources/mirna_hairpin.fasta > resources/mirna.fasta
get_mirna_index : 使用bowtie2
建立miRNA数据的index
bowtie2-build --threads 1 --quiet resources/mirna.fasta resources/mirna > resources/logs/get_mirna_index.log
3.2 下机fastq数据质控
merge_fastq_files : 合并单个样本所有的fastq文件到一个文件中
cat fastq/sampleA_S1_L001_R1_001.fastq.gz \fastq/sampleA_S1_L002_R1_001.fastq.gz \fastq/sampleA_S1_L003_R1_001.fastq.gz \fastq/sampleA_S1_L004_R1_001.fastq.gz > results/sampleA/fastq/sampleA_untrimmed.fastq.gz
trim_illumina_adapters : 从reads的3'端去除Illumina的通用接头序列(Adapter sequences)
cutadapt -a AGATCGGAAGAGCACACGTCTGAAC \--discard-untrimmed \-m 20 -q 10 -j 1 \-o results/sampleA/fastq/sampleA_illumina_trimmed.fastq.gz \--quiet results/sampleA/fastq/sampleA_untrimmed.fastq.gz > results/sampleA/logs/trim_illumina_adapters.log
trim_umi_sequences : 从reads的3'端去除UMI
(Unique Molecular Identifier)序列,并将其添加到reads的头部
umi_tools extract --3prime \--stdin=results/sampleA/fastq/sampleA_illumina_trimmed.fastq.gz \--bc-pattern=NNNNNNNNNNNN \--stdout=results/sampleA/fastq/sampleA_umi_trimmed.fastq.gz > results/sampleA/logs/trim_umi_sequences.log
拓展阅读:单细胞数据分析中的UMI序列
trim_5_adapter : 从reads 5'端去除接头
cutadapt -g GTTCAGAGTTCTACAGTCCGACGATC \-m 20 -j 1 \-o results/sampleA/fastq/sampleA_5adapter_trimmed.fastq.gz \--quiet results/sampleA/fastq/sampleA_umi_trimmed.fastq.gz > results/sampleA/logs/trim_5_adapter.log
trim_3_adapter: 从reads 3'端去除接头,丢弃未修剪的reads并执行一些质量调整
cutadapt -a AACTGTAGGCACCATCAAT \--discard-untrimmed \-m 20 -q 10 -j 1 \-o results/sampleA/fastq/sampleA_trimmed.fastq.gz \--quiet results/sampleA/fastq/sampleA_5adapter_trimmed.fastq.gz > results/sampleA/logs/trim_3_adapter.lo
3.3 比对reads到数据库
map_mirna_sequences: 将reads比对到mirBase数据库序列,直接转换为排序的.bam
文件,并生成.bai
的索引文件。没有比对上的reads将保存为新的fastq
文件,随后与piRNA的数据库比对。
bowtie2 --very-sensitive-local \-p 1 -x resources/mirna \--un-gz results/sampleA/fastq/sampleA_mirna_unmapped.fastq.gz \results/sampleA/fastq/sampleA_trimmed.fastq.gz | samtools sort - > results/sampleA/bam/sampleA_mirna.bam 2> results/sampleA/logs/map_mirna_sequences.log samtools index results/sampleA/bam/sampleA_mirna.bam results/sampleA/bam/sampleA_mirna.bai
map_pirna_sequences: 将reads比对到piRNAdb数据库序列,直接转换为排序的.bam
文件,并生成.bai
的索引文件。有比对上的reads将保存为新的fastq
文件,随后与参考基因组比对(本文没有执行该步骤)。
bowtie2 --very-sensitive-local -p 1 \-x resources/pirna \--un-gz results/sampleA/fastq/sampleA_unmapped.fastq.gz \results/sampleA/fastq/sampleA_mirna_unmapped.fastq.gz | samtools sort - > results/sampleA/bam/sampleA_pirna.bam 2> results/sampleA/logs/map_pirna_sequences.log samtools index results/sampleA/bam/sampleA_pirna.bam results/sampleA/bam/sampleA_pirna.bai
3.4 合并bam文件
merge_bam_files: 将单个 样本BAM 文件合并为一个BAM文件,以减少处理文件的数量并简化后续工作流。
samtools merge --threads 1 \results/sampleA/bam/sampleA_merged.bam \results/sampleA/bam/sampleA_mirna.bam \results/sampleA/bam/sampleA_pirna.bam > results/sampleA/logs/merge_bam_files.log samtools index results/sampleA/bam/sampleA_merged.bam results/sampleA/bam/sampleA_merged.bai
3.5 去除PCR重复
deduplicate_reads: 使用UMItools dedup命令根据先前识别的UMI序列,去除重复的reads
umi_tools dedup \-I results/sampleA/bam/sampleA_merged.bam \-S results/sampleA/bam/sampleA_merged_dedup.bam > results/sampleA/logs/deduplicate_reads.log samtools index results/sampleA/bam/sampleA_merged_dedup.bam results/sampleA/bam/sampleA_merged_dedup.bai
3.6 获得count表
count_deduplicated_read_counts: 统计BAM文件中去除重复reads的计数,并从输出中提取相关列。
samtools idxstats results/sampleA/bam/sampleA_merged_dedup.bam | \awk 'BEGIN{OFS='\t'} $3 != 0' - | \cut -f1,3 -| sed '1 i\Molecule Count' - > results/sampleA/sampleA_counts_dedup.tsv
miRNA生物信息数据分析流程初探(一) - 生信小木屋
snakemake-qiaseq-mirna
↩
Cutadapt Removes Adapter Sequences From High-Throughput Sequencing Reads
↩
UMI-tools: Modelling sequencing errors in Unique Molecular Identifiers to improve quantification accuracy
↩
miRNA生物信息数据分析流程初探(一)相关推荐
- 扩增子和宏基因组数据分析流程和可视化方案—刘永鑫(南京,2020年11月27日)
感谢中科院南京土壤所褚海燕老师的邀请,参加微生物生态与生物信息技术培训. 本次会议预计300人规模的会议,结果现场来了超千人.即使会议进行至第二天下午接近尾声,依然火爆如下: 我将本次90分钟报告&l ...
- 扩增子和宏基因组数据分析流程和可视化方案—刘永鑫(南京,2020年10月27日)
生物信息学习的正确姿势 NGS系列文章包括NGS基础.高颜值在线绘图和分析.转录组分析 (Nature重磅综述|关于RNA-seq你想知道的全在这).ChIP-seq分析 (ChIP-seq基本分析流 ...
- SAP MM 工序委外流程初探
SAP MM 工序委外流程初探 惭愧的说,做D项目之前,笔者都没有玩过SAP里的工序委外功能.因笔者之前参与的项目里,都没有工序委外流程.不过幸运的是,笔者现在所在的D项目里,业务上有工序委外的场景, ...
- producer send源码_Kafka源码深度剖析系列(七)——Producer核心流程初探
本次内容我们有两个目标:第一个初探Producer发送消息的流程第二个我们学习一下Kafka是如何构造异常体系的一.代码分析Producer核心流程初探 //因为生产中开发使用的是异步的方式发送的消息 ...
- QIIME 2:可重复、交互和扩展的微生物组数据分析流程
文章目录 QIIME2:可重复.可交互.适用范围广和可扩展的微生物组数据科学 摘要 正文 图1. 交互式可视化工具 图2. 迭代记录数据来源确保分析可重复 代码可用 在线方法 提取QIIME2的存档内 ...
- python数据分析及可视化(九)pandas数据规整(分组聚合、数据透视表、时间序列、数据分析流程)
作业 拼接多个csv文件 去除重复数据,重新索引 自动挡和手动挡数目 计算每个城市二手车数量 统计每个汽车品牌平均售价价格(不是原价) 分组与聚合 如下表所示,5行3列的表格,5种水果分别对应的名称, ...
- python与数据分析结合_将Python和R整合进一个数据分析流程
原标题:将Python和R整合进一个数据分析流程 在Python中调用R或在R中调用Python,为什么是"和"而不是"或"? 在互联网中,关于"R ...
- 实际业务中的数据分析流程和痛点
平常我们在学校里完成一个数据分析,或者数据挖掘的项目,很多时候的流程是: 在这种分析场景中,我们会更关注如何选择合适的方法来达到我们分析的目的.比如我们现在面对的是一个信用卡欺诈的识别问题,我们已经有 ...
- DPABI详细使用教材——数据准备、预处理流程、数据分析流程
dpabi必看内容 1.DPABI(用于脑成像的数据处理和分析的工具箱)的下载和安装步骤 2.DPABI详细使用教材--数据准备.预处理流程.数据分析流程 3.DPABISurf的安装及使用(wind ...
- 如何做数据分析,数据分析流程是什么?
前言 如何做数据分析,数据分析流程是什么?数据分析是基于商业目的,有目的地进行收集.整理.加工和分析数据,提炼出有价值的信息的一个过程.整个过程大致可分为五个阶段,具体如下图所示. 关于图中流程的相关 ...
最新文章
- Yolo模型部署的两种方法
- Nature:“巨型原子”使芯片同时处理和收发量子信息成为可能
- 制造型企业如何降低成本提升核心竞争力
- SQL改變字符串標識符
- 卷积神经网络 卷积的概念
- python十种日期格式_Python 日期格式相关
- java disjoint_java – Union Find算法的应用(Disjoint Set)
- SQL Server安装问题程序被挂起的错误解决办法
- mysql date max_mysql – 每个ID的SELECT MAX DATE
- python语言的主网址-怎么用Python提取域名中的主域名
- 【渝粤题库】陕西师范大学202281 中央银行学II 作业(专升本)
- 目标高盛,卡方科技用智能金融科技撬动国内量化投资
- 真机测试无法验证应用
- python斐波那契螺旋线怎么画向日葵心_斐波那契螺旋线的图形作法
- 晨枫U盘起动盘制做过程
- 你知道什么是 短路与 和 短路非吗 ???
- INSERT FIRST和INSERT ALL
- Ireport安装使用问题汇总
- pytorch学习资源汇总
- rancher配置日志集中管理
热门文章
- Taro 3.x 开发 APP 记录 (持续记录中。。。)
- jsp 展示服务器pdf文件,jsp实现pdf在线预览功能
- Python函数语法里的中括号和逗号是什么意思
- 微信表情符号写入案件判决
- 分享一本好书《极简主义》-----书中的价值观绝对能引领我们
- SK创新在2019年下半年将实现柔性显示器核心材料FCW量产
- SparkStreaming的背压机制
- linux搭建云存储服务,CentOS 6.3搭建个人私有云存储ownCloud
- 2203-python 24点游戏
- java-php-python-ssm基于智能选课系统的设计与实现计算机毕业设计