python lncrna_LncRNA生信分析案例

LncRNA的筛选

做个笔记，哈哈哈

这个流程是处于链特异性建库(dUTP)的RNA-seq流程，对于非链特异性建库需要修改一些参数。

一、qc和reads质量过滤，fastqc，multiqc，trim_galore就可以了

1.fastqc得到质量报告

ls *gz |xargs -I {} echo 'nohup fastqc {} &'>fastqc.sh

bash fastqc.sh

上面代码比较适合于批量处理，如果电脑核比较紧张，还是用下面这个吧

ls *.fastq.gz | xargs fastqc

样本量大的时候用multiqc，生成汇总信息

multiqc .

2.使用trim_galore进行质量过滤，质量过滤软件很多，随便选一种就好了。

ls *1.fastq.gz > 1

ls *2.fastq.gz > 2

paste 1 2 > config

rm 1 2

bin_trim_galore=trim_galore

dir='/home/dklly/RNA_seq/clean'

cat config |while read pair

do

arr=($pair)

fq1=${arr[0]}

fq2=${arr[1]}

$bin_trim_galore -q 25 --phred33 --length 80 --stringency 3 --paired -o $dir $fq1 $fq2

done

二、比对(hisat2)

1.参考基因组(fasta文件)，注释文件(gtf文件)建议去ensembl下载，ncbi，ucsc上也有，随性就好

2.使用hisat2内置的两个python脚本重gtf文件中获取外显子数据，和剪接位点数据，也可以提取snp信息(需要vcf文件)

hisat2_extract_exons.py

hisat2_extract_splice_sites.py

hisat2_extract_snps_haplotypes_VCF.py

3.建立索引，这个耗时有点长，但是hisat2-build可以使用多线程，呵呵

nohup ~/dk/Miniconda3/envs/RNA_seq_software/bin/hisat2-build -p 8\ Oryctolagus_cuniculus.OryCun2.0.dna_sm.toplevel.fa \

--ss Oryctolagus_cuniculus.OryCun2.0.93.ss\

--exon Oryctolagus_cuniculus.OryCun2.0.93.exon \

Oryctolagus_cuniculus.OryCun2.0.dna_sm.toplevel\

4.对双端reads进行比对，得到sam文件，我直接用samtools转为bam了，节约空间

$hista2 -p 20 --dta --rna-strandness RF -x $index -1 $fq1 -2 $fq2 -S $samName

samtools view -bS -@ 10 $samName > ${samName%.*}.bam

三、转录组的组装

ls *.sort.bam|while read id;do ~/miniconda3/envs/RNA-seq-software/bin/stringtie ${id} -p 2 -G ~/index/Oryctolagus_cuniculus.OryCun2.0.93.gtf -o ${id%%.*}.gtf -l ${id%%.*} ; done

然后使用stringtie --merge合并得到，合并的gft文件

未完待续

收藏 | 0