可变剪切画图
- 创建一个python=2.7的环境
conda create -n rmats402 python=2.7
- 进入环境
conda activate rmats402
- 下载作图软件
conda install -y rmats2sashimiplot
- 根据rmats结果提取想要的基因
- 第一种方法:先挑出想要的基因,记录每个基因的行号和来自哪个文件,使用sed命令提取这些行号
sed -n '1p;552p;553p;984p;1382p;1980p;2590p;4064p;4655p;4883p;4884p;6222p;6500p;6621p;7534p;7870p;7871p;7872p;8953p;8962p;11164p;11989p;12544p;13226p;13474p;13577p;14115p;14116p;15476p;17042p;17043p;17517p' SE.MATS.JC.txt > ../../YG.SE.MATS.txt
- 第二种方法: 记录每个基因ID号,使用awk命令提取带有ID号的这些行
cat A3SS.MATS.JC.txt | awk '$1=="ID" || $1==13 {print$0}' > ../../YG.A3SS.MATS.txt
cat A5SS.MATS.JC.txt | awk '$1=="ID" || $1==291 {print$0}' > ../../YG.A5SS.MATS.txt
cat MXE.MATS.JC.txt | awk '$1=="ID" || $1==1326 || $1==1484 || $1==1490 {print$0}' > ../../YG.MXE.MATS.txt
cat SE.MATS.JC.txt | awk '$1=="ID" || $1==550 || $1==551 || $1==982 || $1==1380 || $1==1978 || $1==2588 || $1==4062 || $1==4653 || $1==4881 || $1==4882 || $1==6220 || $1==6498 || $1==6619 || $1==7532 || $1==7868 || $1==7869 || $1==7870 || $1==8951 || $1==8960 || $1==11162 || $1==11987 || $1==12542 || $1==13224 || $1==13472 || $1==13575 || $1==14113 || $1==14114 || $1==15474 || $1==17040 || $1==17041 || $1==17515 {print$0}' > ../../YG.SE.MATS.txt
- 处理sort.bam文件 不处理也会报错
samtools view -H H1_sort.bam | sed -e 's/SN:\([0-9XY]\)/SN:chr\1/' -e 's/SN:MT/SN:chrM/' | samtools reheader - H1_sort.bam > H1_sort_chr.bam
samtools view -H H2_sort.bam | sed -e 's/SN:\([0-9XY]\)/SN:chr\1/' -e 's/SN:MT/SN:chrM/' | samtools reheader - H2_sort.bam > H2_sort_chr.bam
samtools view -H H4_sort.bam | sed -e 's/SN:\([0-9XY]\)/SN:chr\1/' -e 's/SN:MT/SN:chrM/' | samtools reheader - H4_sort.bam > H4_sort_chr.bam
samtools view -H L1_sort.bam | sed -e 's/SN:\([0-9XY]\)/SN:chr\1/' -e 's/SN:MT/SN:chrM/' | samtools reheader - L1_sort.bam > L1_sort_chr.bam
samtools view -H L3_sort.bam | sed -e 's/SN:\([0-9XY]\)/SN:chr\1/' -e 's/SN:MT/SN:chrM/' | samtools reheader - L3_sort.bam > L3_sort_chr.bam
samtools view -H L2_sort.bam | sed -e 's/SN:\([0-9XY]\)/SN:chr\1/' -e 's/SN:MT/SN:chrM/' | samtools reheader - L2_sort.bam > L2_sort_chr.bam
- 运行命令SE格式
rmats2sashimiplot --b1 H1_sort_chr.bam,H2_sort_chr.bam,H4_sort_chr.bam --b2 L1_sort_chr.bam,L2_sort_chr.bam,L3_sort_chr.bam -t SE -e ./YG.SE.MATS.txt --l1 H --l2 L --exon_s 1 --intron_s 1 -o YG_SE_MATS_plot
- 运行命令 A5SS格式
rmats2sashimiplot --b1 H1_sort_chr.bam,H2_sort_chr.bam,H4_sort_chr.bam --b2 L1_sort_chr.bam,L2_sort_chr.bam,L3_sort_chr.bam -t A5SS -e ./YG.A5SS.MATS.txt --l1 H --l2 L --exon_s 1 --intron_s 1 -o YG_A5SS_MATS__plot
- 运行命令 A3SS格式
rmats2sashimiplot --b1 H1_sort_chr.bam,H2_sort_chr.bam,H4_sort_chr.bam --b2 L1_sort_chr.bam,L2_sort_chr.bam,L3_sort_chr.bam -t A3SS -e ./YG.A3SS.MATS.txt --l1 H --l2 L --exon_s 1 --intron_s 1 -o YG_A3SS_MATS__plot
- 成功运行
Reading sample label: PLD5 L-2 IncLevel: 0.00
Skipping read with CIGAR 38M1I111M
Skipping read with CIGAR 36M1I113M
Skipping read with CIGAR 36M1I113M
Skipping read with CIGAR 112M2D38M
Skipping read with CIGAR 49M2D101M
Skipping read with CIGAR 130M1D20M
Skipping read with CIGAR 23M1D127M
Reading sample label: PLD5 L-3 IncLevel: 0.00
Skipping read with CIGAR 11M1I138M
Skipping read with CIGAR 30M1I119M
Skipping read with CIGAR 135M2D15M
Skipping read with CIGAR 94M1I55M
Skipping read with CIGAR 96M2D54M
Skipping read with CIGAR 51M2D99M
Skipping read with CIGAR 113M1I36M
Skipping read with CIGAR 113M1I36M
Skipping read with CIGAR 113M1I36M
Skipping read with CIGAR 1S112M1I36M
Saving plot to: /public/jychu/YangGe/cleandata/bam/YG_SE_MATS_plot/Sashimi_plot/chr3:35573451:35573636:+@chr3:35606603:35606681:+@chr3:35637366:35637502:+.pdf
- 如果是自己创建输入文件,则会出现下面的bug
File "/home/ #/anaconda2/lib/python2.7/site-packages/MISO/misopy/sashimi_plot/sashimi_plot.py", line 293, in <module>
main()
File "/home/#/anaconda2/lib/python2.7/site-packages/MISO/misopy/sashimi_plot/sashimi_plot.py", line 289, in main
plot_label=plot_label)
File "/home/#/anaconda2/lib/python2.7/site-packages/MISO/misopy/sashimi_plot/sashimi_plot.py", line 148, in plot_event
%(event_name, pickle_dir)
.....
切记不要用EXCEL创建输入文件
网友评论