一.数据处理
1.得到exon.pos;3UTR.pos;5UTR.pos;repeat.pos文件
##进入小环境运行
conda create -n RNA_editing
conda activate RNA_editing
# 得到exon.pos
grep "exon" Iyun.gff > Iyun.exon.gff
awk -F ';' '{print $1}' Iyun.exon.gff > Iyun.exon.gff1
awk '{print $9"\t"$1"\t"$7"\t"$4"\t"$5}' Iyun.exon.gff1 > Iyun.exon.gff2
sed 's/.exon/_CDS-N/' Iyun.exon.gff2 > exon.pos
sed 's/ID=//g' exon.pos > exon1.pos
#得到repeat.pos
sed 's/Hic_asm_0/Chr01/g' All_Repeat.gff | sed 's/Hic_asm_10/Chr11/g' | sed 's/Hic_asm_1/Chr02/g' | sed 's/Hic_asm_2/Chr03/g' | sed 's/Hic_asm_3/Chr04/g' |sed 's/Hic_asm_4/Chr05/g' | sed 's/Hic_asm_5/Chr06/g' |sed 's/Hic_asm_6/Chr07/g' |sed 's/Hic_asm_7/Chr08/g' |sed 's/Hic_asm_8/Chr09/g' | sed 's/Hic_asm_9/Chr10/g' > All_Repeat1.gff
awk -F ';' '{print $1}' All_Repeat1.gff > All_Repeat1.gff1
awk '{print $9"\t"$1"\t"$7"\t"$4"\t"$5}' All_Repeat1.gff1 > repeat.pos
#得到UTR.pos
grep "three_prime_UTR" Iyun.gff > Iyun.3UTR.gff
awk -F ';' '{print $1}' Iyun.3UTR.gff > Iyun.3UTR.gff1
awk '{print $9"\t"$1"\t"$7"\t"$4"\t"$5}' Iyun.3UTR.gff1 > Iyun.3UTR.gff2
sed 's/.utr3p/_UTR3-N/' Iyun.3UTR.gff2 > Iyun.3UTR.gff3
sed 's/ID=//g' Iyun.3UTR.gff3 > Iyun.3UTR.gff4
#得到ncRNA文件
cat miRNA.gff3 rRNA.gff3 snRNA.gff3 tRNA.gff3 > ncRNA.gff3
awk -F ';' '{print $1}' ncRNA.gff3 > ncRNA1.gff3
awk '{print $9"\t"$1"\t"$7"\t"$4"\t"$5}' ncRNA1.gff3 > ncRNA2.gff3
sed 's/ID=//g' ncRNA2.gff3 > ncRNA.pos
2.得到intron.pos文件
根据上述gtf文件进行后续处理
2.1.得到intron.bed文件
1,从gtf中提取外显子
awk 'BEGIN{OFS="\t";} $3=="exon" {print $1,$4-1,$5}' file.gtf |bedtools sort |bedtools merge -i - > my_exon.bed
2. 提取内含子
awk 'BEGIN{OFS="\t";} $3=="transcript" {print $1,$4-1,$5}' file.gtf |bedtools sort |bedtools subtract -a stdin -b my_exon.bed > my_intron.bed
2.2.对intron.bed进行处理,得到含有geneid的intron.bed文件
此时得到的的内含子文件只有三列,需要把这三列和geneID联系起来
如图所示,my_intron.bed的第二列,和 exon.pos的第五列是对应的,所以可以从exon.pos中提取my_intron.bed相对应的geneid.
image.png
image.png
2.2.1 处理intron.bed文件,得到my_intron_result.txt
#.为防止不同染色体上有相同位置,把my_intron.bed第二列加上对应的Chr,如图
awk '{print $1}' my_intron.bed > 1.txt
awk '{print $2}' my_intron.bed > 2.txt
paste -d '_' 1.txt 2.txt >my_intron_result.txt
image.png
2.2.2处理exon.pos文件,得到CDS_chr.txt
awk '{print $2}' exon.pos > 3.txt
awk '{print $5}' exon.pos > 4.txt
paste -d '_' 3.txt 4.txt > 5.txt
awk '{print $1,$2,$3}' exon.pos > 6.txt
paste 6.txt 5.txt > CDS_chr.txt
rm 1.txt 2.txt 3.txt 4.txt 5.txt 6.txt
image.png
2.2.3 批量从CDS_chr.txt中提取my_intron_result.txt对应的geneid
awk '{print $1}' /home/lx_sky6/yyj/RNA_editing/Iyun/data/intron_all/my_intron_result.txt |xargs -i echo awk "'\$4==\"{}\" {print \$1,\$2,\$3,\$4}'" /home/lx_sky6/yyj/RNA_editing/Iyun/data/intron_all/CDS_chr.txt ">" /home/lx_sky6/yyj/RNA_editing/Iyun/data/intron_all/1.result/{}_intron.bed > 1.pbs
image.png
随后,把提取出来的文件放在一起,检查一下是否与my_intron.bed的数量相等
cat Chr01_* Chr02_* Chr03_* Chr04_* Chr05_* > 1.5.bed
cat Chr06_* Chr07_* Chr08_* Chr09_* Chr10_* Chr11_* > 6.11.bed
cat 1.5.bed 6.11.bed scaffold* > all.bed
cat all.bed |wc -l
cat my_intron.bed | wc -l
如果相等,无误后直接用就可
如果不相等,需要把my_intron.bed文件中重复的进行剔除
awk '{print $4}' all.bed > all_4.txt
sort all_4.txt | uniq -d
把上述找到的重复的,在all.bed里进行剔除(手工),再次检查 all.bed与my_intron_uniq.bed数量是否相同。相同即可把两个文件paste在一起。
paste all.bed my_intron_uniq.bed > intron1.bed
检查all.bed my_intron_uniq.bed合并在一起的文件是否相对应,对应即进行如下操作。
awk '{print $1,$2,$3,$6,$7}' instron3.bed > intron_try.pos
sed 's/_CDS-N/_intron-N/g' intron.pos > intron.pos
网友评论