可变剪切

作者: 宗肃書 | 来源:发表于2021-01-21 17:56 被阅读0次

    可变剪切画图

    • 创建一个python=2.7的环境
    conda create -n rmats402 python=2.7
    
    • 进入环境
    conda activate rmats402
    
    • 下载作图软件
    conda install -y rmats2sashimiplot
    
    • 根据rmats结果提取想要的基因
    1. 第一种方法:先挑出想要的基因,记录每个基因的行号和来自哪个文件,使用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
    
    1. 第二种方法: 记录每个基因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创建输入文件

    相关文章

      网友评论

        本文标题:可变剪切

        本文链接:https://www.haomeiwen.com/subject/yxydzktx.html