Ribo-seq的下游分析方法1-ORFquant以及RiboQC

学习文件：https://github.com/lcalviell/ORFquant，数据也是这篇文章里面的。
1.安装软件

library(devtools)
install_github(repo = "lcalviell/ORFquant")
library(ORFquant)
library(RiboseQC)
conda install -c conda-forge pandoc

安装过程由于是在linux里面安装的R，所以有很多附属的R包没有安装进去，可以用bioconductor也可以用conda进行安装。

2.prepare_annotation_files

?prepare_annotation_files

在R中运行，会出现所需要准备的文件并且会有相应的要求，例如：

因此
根据要求，写进去自己的数据
(1)设置文件夹

getwd（）##获取现在所在的位置
dir.create("1.my_ORFquant_dir")###创建一个新的文件夹
setwd("1.my_ORFquant_dir")###进入这个新创建的文件夹

（2）把fasta转换成2bit格式
用conda下载这三个软件faToTwoBit, twoBitInfo, and twoBitToFa
在genecode里面下载GTF和FASTA文件

faToTwoBit GRCh38.p13.genome.fa GRCh38.p13.genome.2bit

（3）赋值文件

my_riboseq_bam <- "./SRR2433794.bam"
my_gtf_file <- "./gencode.v40.annotation.gtf"
my_fasta_file <- "./GRCh38.primary_assembly.genome.2bit"
####查看是否负值成功
stopifnot(
file.exists(my_riboseq_bam)&
file.exists(my_gtf_file)&
file.exists(my_fasta_file)
)

（4）生成annotation_files

orfquant_anno_file <-prepare_annotation_files(annotation_directory = "annotation_directory/",twobit_file = my_fasta_file,gtf_file = my_gtf_file,scientific_name = "Homo.sapiens",annotation_name = "genc25",export_bed_tables_TxDb = T,forge_BSgenome = Tcreate_TxDb = T)

这一步结束会得到很多文件，大小如下所示

drwxr-xr-x 5 med-zhouh med-chenh 4.0K Jul  4 19:58 BSgenome.Homo.sapiens.genc25
-rw-r--r-- 1 med-zhouh med-chenh 8.3M Jul  4 20:18 cds_txs_coords_similbed.bed
-rw-r--r-- 1 med-zhouh med-chenh 2.4M Jul  4 20:18 fiveutrs_similbed.bed
-rw-r--r-- 1 med-zhouh med-chenh  61M Jul  4 20:18 gencode.v40.annotation.gtf_Rannot
-rw-r--r-- 1 med-zhouh med-chenh 127M Jul  4 20:04 gencode.v40.annotation.gtf_TxDb
-rw-r--r-- 1 med-zhouh med-chenh 4.7K Jul  4 20:18 genetic_codes
-rw-r--r-- 1 med-zhouh med-chenh  490 Jul  4 19:57 GRCh38.primary_assembly.genome.2bit_Homo.sapiens_seed
-rw-r--r-- 1 med-zhouh med-chenh  15M Jul  4 20:18 introns_similbed.bed
-rw-r--r-- 1 med-zhouh med-chenh 3.4M Jul  4 20:18 ncIsof_similbed.bed
-rw-r--r-- 1 med-zhouh med-chenh 4.7M Jul  4 20:18 ncRNAs_similbed.bed
-rw-r--r-- 1 med-zhouh med-chenh 4.3K Jul  4 20:18 seqinfo
-rw-r--r-- 1 med-zhouh med-chenh  18M Jul  4 20:18 table_gene_tx_IDs
-rw-r--r-- 1 med-zhouh med-chenh 2.0M Jul  4 20:18 threeutrs_similbed.bed

（5）QC

#Run the RiboseQC pipeline
ribseqc_file <- RiboseQC_analysis(annotation_file = orfquant_anno_file,bam_files = my_riboseq_bam,fast_mode = T,create_report = T)
#check that the file of p-site offsets has been created
my_offset_file <- paste0(my_riboseq_bam,"_P_sites_calcs")
rl_cutoff_file <- read.delim(my_offset_file)
rl_cutoff_file <- rl_cutoff_file[,c(1,3,9)]
write.table(rl_cutoff_file,file = 'rl_cutoff_file',quote = FALSE,sep = '\t',row.names = FALSE)
#verify our psite offset file exists
file.exists(my_offset_file)

生成很多文件

 SRR2433794.bam_P_sites_minus.bedgraph
SRR2433794.bam_P_sites_plus.bedgraphSRR2433794.bam_P_sites_uniq_minus.bedgraph
SRR2433794.bam_P_sites_uniq_plus.bedgraph
SRR2433794.bam_coverage_minus.bedgraph       SRR2433794.bam_results_RiboseQC
SRR2433794.bam_coverage_plus.bedgraph        SRR2433794.bam_results_RiboseQC_all
SRR2433794.bam_coverage_uniq_minus.bedgraph  SRR2433794.bam_RiboseQC_report_files
SRR2433794.bam_coverage_uniq_plus.bedgraph   SRR2433794.bam_RiboseQC_report.html_plots
SRR2433794.bam_for_ORFquant                  SRR2433794.bam_RiboseQC_report.html_report_text_output.txt
SRR2433794.bam_for_SaTAnn
SRR2433794.bam_junctions
SRR2433794.bam_P_sites_calcs

（6）准备P-site文件和运行ORFquant

my_orfquant_psites_file < prepare_for_ORFquant(orfquant_anno_file,my_riboseq_bam,path_to_rl_cutoff_file = 'rl_cutoff_file')
#run orfquant
ORFquant_results <- run_ORFquant(for_ORFquant_file = my_orfquant_psites_file,annotation_file = orfquant_anno_file,n_cores = 40)
####在超算里面运行核数n可以用到40多，但是如果是本地就只能1个核，在本地时间大概12h左右##
orfquant_res_file <- paste0(my_riboseq_bam,"_for_ORFquant_final_ORFquant_results")
stopifnot(file.exists(orfquant_res_file))

(7)plot ORFquant results

plot_ORFquant_results(for_ORFquant_file=my_orfquant_psites_file,ORFquant_output_file=orfquant_res_file, annotation_file=orfquant_anno_file)
plotfolder <- paste0 (orfquant_res_file,"_plots/")
stopifnot (file.exists(plotfolder))
orfquantfile=paste0("Ribo-WT-23Aligned.sortedByCoord.out.","bam_for_ORFquant_final_ORFquant_results_plots/","Ribo-WT-23Aligned.sortedByCoord.out.bam_for_ORFquant_ORFquant_plots_RData")
create_ORFquant_html_report(input_files = orfquantfile,input_sample_names = "Ribo-WT-23Aligned.sortedByCoord.out.",output_file= "Ribo-WT-23Aligned.sortedByCoord.out_ORFquant_report.html")- GTF和FASTA一定要是对应的。- 比对过程不要改动参数。- bam文件推荐用star进行比对。- 祝成功！

Ribo-seq的下游分析方法1-ORFquant以及RiboQC相关推荐

从测序到宏基因组：聚焦菌群生信分析方法最前沿
今天是第2439期日报. Nature子刊:基于三代测序的宏基因组分析助力完整微生物进化研究 Nature Methods[IF:47.99] ① 基于Pacific Biosciences或Oxfo ...
结合语言知识和深度学习的中文文本情感分析方法
结合语言知识和深度学习的中文文本情感分析方法徐康庭, 宋威北方工业大学信息学院摘要:在目前的中文文本情感分析研究中,基于语义规则和情感词典的方法通常需要人工设置情感阈值:而基于深度学习的方法由于 ...
易基因｜综述：单细胞DNA甲基化分析方法全介绍及未来发展前景预测
大家好,这是专注表观组学十余年,领跑多组学科研服务的易基因. 2021年07月10日,<Biomolecules>杂志上发表一篇关于单细胞表观测序的综述文章,详细介绍了单细胞DNA甲基化 ...
转录组解读及下游分析
普通真核转录组解读约40min , 原核转录组的大部分分析点与普通真核转录组一致. 易基因 - 转录组结题报告讲解_哔哩哔哩_bilibili 分析软件准备需要准备好 Tbtools 软件,后续大 ...
易基因：全基因组CpG密度和DNA甲基化分析方法比较（MeDIP、RRBS和WGBS）| 研究综述
大家好,这里是专注表观组学十余年,领跑多组学科研服务的易基因. CpG密度(CpG density)与各种组织中的DNA甲基化相关.基因组按CpG密度分为:CpG岛(CpG island,CGI).C ...
微生物差异丰度分析方法的综合评价：现状及潜在解决方案
微生物差异丰度分析方法的综合评价:现状及潜在解决方案原文〉A comprehensive evaluation of microbial differential abundance analysi ...
面向过程（或者叫结构化）分析方法与面向对象分析方法到底区别在哪里？
AutoSAR入门到精通系列讲解将从2019年开始更新关于AutoSAR的知识,从入门到精通,博主xyfx和大家一起进步雪云飞星 ¥29.90 去订阅简单地说结构化分析方法主要用来分析系统的功能 ...
java三维滑雪,第六章三维数据空间分析方法.ppt
第六章三维数据空间分析方法 * * * * 可视性分析 * * 自然邻域法插值基本思路: 利用输入点及邻近栅格单元进行插值生成栅格表面. 方法: 利用输入数据点(样本点)为节点,建立Delauna ...
mysql慢查询开启及分析方法
最近服务维护的公司的DB服务器,总是会出现问题,感觉需要优化一下了,登陆上去,发现慢查询日志都没有开,真是惭愧, 故果断加上慢查询日志,经过分析sql记录,发现问题很多,开发人员很多没有对sql优化, ...

Ribo-seq的下游分析方法1-ORFquant以及RiboQC

Ribo-seq的下游分析方法1-ORFquant以及RiboQC相关推荐

最新文章

热门文章