鉴定新的lncRNA之上游流程
上游流程，通常指的是ngs测序数据fastq文件，在服务器级别的计算资源里面的一系列处理。因为个人电脑很难hold住，而且流程很少变动，所以通常是公司代替客户完成。属于吃力不讨好的技能，学习成本极高，但是学完之后用的频率很低！下游分析流程，指的是拿到上游分析结果之后在自己个人电脑里面的统计可视化图表制作。
比如RNA-seq数据，上游就是fastq的质量控制，比对，定量，最后拿到表达矩阵。而下游就是表达矩阵的一系列统计学分析， 包括PCA,相关性热图，层次聚类图，差异分析，火山图，表达量热图，GO/KEGG数据库功能注释等等。
对lncRNA-seq鉴定新的lncRNA流程来说，也是上下游独立。我们先介绍一下上游流程哦：
首先使用conda来创建LncRNA-seq的实战软件环境
conda create -n lncRNA
conda activate lncRNA
conda install -y -c bioconda hisat2 stringtie samtools fastp gffcompare
# conda search gffcompare然后下载猪的参考基因组fasta序列，并且构建hisat2的索引文件
mkdir -p ~/reference/genome/
cd ~/reference/genome/
mkdir pig
cd pig
wget ftp://ftp.ensembl.org/pub/release-99/fasta/sus_scrofa/dna/Sus_scrofa.Sscrofa11.1.dna.toplevel.fa.gz
gunzip -d Sus_scrofa.Sscrofa11.1.dna.toplevel.fa.gz
# 需要注意文件大小，以及参考基因组是否下载成功哦！
wget ftp://ftp.ensembl.org/pub/release-99/gtf/sus_scrofa/Sus_scrofa.Sscrofa11.1.99.chr.gtf.gz
gunzip -d Sus_scrofa.Sscrofa11.1.99.chr.gtf.gz
nohup hisat2-build -p 4 Sus_scrofa.Sscrofa11.1.dna.toplevel.fa pig_hisat2 & 一定要注意 hisat2-build 根据参考基因组的fasta序列文件构建索引对内存的消耗，以及输出文件的大小哦。
接着下载文章lncRNA-seq的fastq测序数据
cat fq.txt |while read id
do
ascp -QT -l 300m -P33001 \
-i ~/miniconda3/envs/download/etc/asperaweb_id_dsa.openssh \
era-fasp@$id .
done
# nohup bash step1-aspera.sh 1>step1-aspera.log 2>&1 &这个脚本会根据你在EBI里面搜索到的 fq.txt 路径文件，来批量下载fastq测序数据文件。这个 fq.txt 文件里面的路径还是蛮有规律的，如下：
fasp.sra.ebi.ac.uk:/vol1/fastq/SRR180/009/SRR1805929/SRR1805929_1.fastq.gz
fasp.sra.ebi.ac.uk:/vol1/fastq/SRR180/009/SRR1805929/SRR1805929_2.fastq.gz
fasp.sra.ebi.ac.uk:/vol1/fastq/SRR180/000/SRR1805930/SRR1805930_1.fastq.gz
fasp.sra.ebi.ac.uk:/vol1/fastq/SRR180/000/SRR1805930/SRR1805930_2.fastq.gz值得注意的是，EBI数据库也不是那么稳定，有时候这个下载会失败，大家可以等候几天重新尝试，或者干脆转到去sra数据库下载。
批量运行fastp+hisat2+stringtie流程
需要学习和使用3个软件的用法，fastp+hisat2+stringtie，理解参数：
conda activate lncRNA
index=/home/yb77613/reference/genome/pig/pig_hisat2
gtf=/home/yb77613/reference/genome/pig/Sus_scrofa.Sscrofa11.1.99.chr.gtf
fastp -i 1.raw_fq/${id}_1.fastq.gz \
-o 2.clean_fq/${id}_1.fastp.fq.gz \
-I 1.raw_fq/${id}_2.fastq.gz \
-O 2.clean_fq/${id}_2.fastp.fq.gz \
-l 36 -q 20 --compression=6 \
-R ${id} -h ${id}.html
fq1=2.clean_fq/${id}_1.fastp.fq.gz
fq2=2.clean_fq/${id}_2.fastp.fq.gz
hisat2 -p 4 -x $index -1 $fq1 -2 $fq2 | \
samtools sort -@ 4 -o 3.hisat2_bams/$sample.bam -
stringtie -p 4 -G $gtf \
-o 4.stringtie_gtfs/$sample.gtf \
-l $sample 3.hisat2_bams/$sample.bam 批量处理！
合并gtf后鉴定和筛选新的lncRNA
conda activate lncRNA
ls `pwd`/4.stringtie_gtfs/*gtf > 5.lncRNA/gtf.list
cd 5.lncRNA
gtf=/home/yb77613/reference/genome/pig/Sus_scrofa.Sscrofa11.1.99.chr.gtf
nohup stringtie --merge -p 6 -G $gtf -o stringtie_merged.gtf gtf.list > merge.log 2>&1 &
gff=Sus_scrofa.Sscrofa11.1.93.gff3
gffcompare -r $gff -o gffcomp93 stringtie_merged.gtf
# 49448 reference transcripts loaded.
# 69 duplicate reference transcripts discarded.
# 113472 query transfrags loaded.
gtf=/home/yb77613/reference/genome/pig/Sus_scrofa.Sscrofa11.1.99.chr.gtf
gff=$gtf
gffcompare -r $gff -o gffcomp99 stringtie_merged.gtf
# 60954 reference transcripts loaded.
# 84 duplicate reference transcripts discarded.
# 113472 query transfrags loaded.可以看到，使用不同版本的gtf文件，对结果的影响很明显哦，鉴定新的lncRNA的话，也许2016年你认为是新的lncRNA，到2018年就已经被记录在册啦。
如果是人类数据，我们通常是使用gencode的gtf文件哦
gtf=/home/yb77613/reference/gtf/gencode/gencode.v32.annotation.gtf
nohup stringtie --merge -p 6 -G $gtf -o stringtie_merged.gtf gtf.list > merge.log 2>&1 &
gffcompare -r $gtf -o gffcomp stringtie_merged.gtf其它物种，大家类推即可，主要是人类和小鼠，做这个流程的意义不是很大。当然，如果你强行跑流程，也是可以的，结果的生物学解释就比较考验人。
后面的ngs下游流程更精彩：
举报
腾讯云开发者
