2018-09-03

第一天 数据下载
https://sra-explorer.info/?

SRA Explorer

流程：全部选中--add to collection--saved datasets--Aspera commands for downloading FastQ files

第二天 数据过滤
安装conda
流程：下载minconda--配置仓库（condarc）--安装软件
htop F5显示进程
下载fastp
替换”rename 's///' *.gz“
打印文件名ls *.gz >test
打印文件中第一列 awk '{print 1 -o “clean/"1".html -j"$1".json &" }' +file >run_fastp.sh

第三天 参考基因组下载
Ensemble下载基因组序列，注释文件，蛋白序列
基因组序列：下载primary_assembly
toplevel里面包含搞不明白的DNA序列，不考虑
下载的数据，解压 cat *.fa>genome.fa

注释文件：chr.gtf全部定位到基因组了
abinitio从头注释，软件预测，不靠谱
将gff转变为gtf代码：gffread -T -o file.gtf file.gff3
蛋白文件：awk '{print 1}' file
第三列中基因的个数 awk '$3=="gene"' file |wc -l

第四天

软件包FastX-toolkit用于裁断基因每行的个数，80个，长度转化
conda install fastx-toolkit
一行70个字符：fasta_formatter -i genome.fa -o genome_formater.fa -w 70

分析思路：
比对到基因组： Hisat2，STAR
比对到转录组; Bowties

饱和性曲线Saturation curve
6G：碱基数目 对于大多数动植物样本，推荐20M reads X （150+150）=6000M

表达定量Quantification
subread --featureCounts进行定量
链特异性文库可以区分正义量

标准化Normalization
TPM

比对
第一步（构建参考基因组）：参考基因组处理，构建index，排序。用于比对顺利
输入：基因组序列 输出：构建好参考基因组 软件：hisat2
代码：`hisat2-build ../ref/genome.fasta ../ref/genome 1>hisat2-build.log 2>&1
第二步（比对）：
输入：构建好的基因组，测序数据 输出：比对结果（sam）软件：hisat2
代码

hisat2 --new-summary -p 10 -x ../ref/genome -U ../data/BLO_S1_LD1.fq.gz -S BLO_S1_LD1.sam --rna-strandness R 1>BLO_S1_LD1.log 2>&1
第三步（压缩和排序）：
输入：sam 输出：bam 软件：samtools sort
代码
samtools sort -o BLO_S1_LD1.bam BLO_S1_LD1.sam
第四步：bam文件索引
输入：bam 输出bam.bai 软件：
代码:samtools index BLO_S2_LD3.bam

IGV：基因组文件 基因注释文件 bam文件 sample.bam.bai
第一步：构建基因组：输入基因组文件，基因组解释文件。

基因组构建

第二步：将bam文件，bam.bai文件放于一个文件夹
在IGV中选择File--loadfile--选择bam文件即可，无需选择bam.bai文件

第五天 定量

multiqc生成比对报告
samtools flagstat file.bam 查看比对情况
log文件更推荐查看比对

IGV查看比对情况：
HIsat2--stingtie：通过转录组拼接对于原来的基因结构优化（不推荐）
PASA对原来的基因结构进行优化

个体重测序本质是基因分型
群体重测序本质是等位基因频率

表达定量：
输入得到的bam文件，基因注释文件，输出：定量结果文件，软件：Rsubread（featurecounts）R语言
安装：(1)R包来源于于CRAN
install.packages("")
(2)R包来源于Bioconductor，查看安装文档

Rscript file.R
代码：Rscript script/run-featurecounts.R -b ../1.Mapping/BLO_S1_LD1.bam -g ../ref/genes.gtf -o BLO_S1_LD1

-b:bam文件 -g:gtd文件 -o：输出文件前缀

批量脚本：awk '{print "Rscript script/run-featurecounts.R -b ../1.Mapping/" $2".bam -g ../ref/genes.gtf -o "$2}' ../data/samples.txt

生成count

第六天 表达定量生成矩阵