02

保守的 isforms鉴定

使用ab-Blast进行鉴定

使用xdformat建库

  • -n指定序列类型

  • -o 输出数据库前称

  • 最后跟fasta文件

xdformat  -n  -o test test.fa

进行比对

参考说明 https://blast.advbiocomp.com/doc/parameters.html

  • 库前称

  • 序列文件

  • -o输出文件

~/software/wu-blast/ab-blast-20200317-linux-x64/ab-blastn test test.fa -o 11  e=1e-5 mformat=2

提取isform每个剪切位点附近的exon序列

  • 提取exon两端锚点附近300bp序列

  • 不足300bp则提取整个exon

运行wu-blast

  • 同源基因都存在isform的情况

python ~/scripte/Alternative/module/FEST3/extract_spliceSiteSeq.py -p geneExpress/isform.gtf  -r gene_isformCount.txt2  -g ~/work/Alternative/data/Ga_genome/G.arboreum.Chr.v1.0.fasta -o tmp
## 建库
~/software/wu-blast/ab-blast-20200317-linux-x64/xdformat -n -o all_vs_all all_vs_all.fa
##鉴定同源基因间保守的isform

python ~/scripte/Alternative/module/FEST3/wu-Blast.py  -homolog ~/work/Alternative/result/homologo/homologGene/A2_vs_At_collinearity.txt  -r1 ../A2/gene_isformCount.txt2  -r2 ../TM1/At_isformCount.txt  -fa1 A2_isform_splice.fa  -fa2 At_isform_splice.fa  -o hhhhhh
##

鉴定保守isform2.0

  • exon数目相同，

  • 从5‘到3’端每个exon的长度都接近

python ~/scripte/Alternative/module/FEST3/conserveIsform.py  -r1 ../A2/gene_isformCount.txt2  -r2 ../D5/gene_isformCount.txt2 -gtf1 ~/work/Alternative/result/Ga_result/CO11_12_result/07_annotation/merge.gtf -gtf2 ~/work/Alternative/result/Gr_result/CO41_42_result/07_annotation/merge.gtf -gap 15 -homolog ~/work/Alternative/result/homologo/homologGene/A2_vs_D5_collinearity.txt  -o zzz

鉴定保守的isform3.0

  • wu-blast结果存在多对多的情况

  • 比较isform的长度，取最接近的那个

  • 去除isform长度差异不超过200bp的

##转换wu-blast结果
awk '$1~/evm/{print $0}$1~/Ghir/{print $2"\t"$1}' ../A2_At.txt
##提取保守的isform
python ~/scripte/Alternative/module/FEST3/conserveIsform2.py  -blast A2_At.txt  -r1 ~/work/Alternative/result/Ga_result/CO11_12_result/07_annotation/merge.gtf  -r2 ~/work/Alternative/result/Gh_result/CO31_32_result/07_annotation/merge.gtf  -o 222
##去重
sort 222 |uniq |sort -k2,2 -k5,5n|awk '{print $1,$3,$4,$5,$2}' OFS="\t" |uniq -f4|awk '{print $5,$2,$3,$4,$1}' OFS="\t" |sort -k5,5 -k4,4n|uniq -f4 |awk '{print $5,$1,$2,$3,$4}' OFS="\t" >A2_At.txt
## 筛选并且更改格式
awk '$5<=200{split($1,a,"-");split($2,b,"-");print a[1]"\t"a[2]"\t"b[1]"\t"b[2]}' ../conserve_isform/conserve_isform/D5_Dt.txt

保守的isform统计

比较

都存在isform的基因数

isform数

保守isform数

基因数

A2 D5

13543

77180

11823

8196

At Dt

10368

50282

6331

4850

A2 At

12041

68419

12868

8604

D5 Dt

12580

63143

12335

8570

统计存在保守isfrom的数目分布，以及对应的gene isform的总数的平均值

python ~/scripte/Alternative/module/FEST3/isformlevel/isformNumber.py -c conserve_isform.txt  -all isform_count.txt  -f1 0 -f2 2 -o1 test -o2 test2

二倍体到四倍体，多倍化过程中，存在保守isform的基因，非保守的isform减少，表明isform的多态性在减少

python ~/scripte/Alternative/module/FEST3/isformConserveRate.py -all isform_count.txt  -c conserve_isform.txt  -o 11

二倍体内保守，四倍体内不保守例子

isoform在多倍化后的丢失

  • A2 与D5保守，只有At存在保守的

  • A2与D5保守，只有Dt存在保守的

二倍体内部保守，四倍体内保守的例子

isoform多倍化后，在两个亚基因组间共享

二倍体内保守，四倍体内同样保守

python ~/scripte/Alternative/module/FEST3/isformlevel/polypolid.py -homolog ~/work/Alternative/result/homologo/homologGene/A2_D5_At_Dt_collinearity.txt  -A2D5 ../A2_vs_D5/conserve_isform.txt  -AtDt ../At_vs_Dt/conserve_isform.txt  -A2At ../A2_vs_At/conserve_isform.txt  -D5Dt ../D5_vs_Dt/conserve_isform.txt  -o AtDt分化 -o1 AtDt趋同 -o2 保持不变

  • 输出的结果还需要去一下重,由于在用wu-blast直接比对的时候没有鉴定出保守性，但是通过两个基因组间的关系间接证明了保守性，这个先去掉，原始文件先存着；到时候设计引物的时候在看看

## 多倍化后isoform丢失
awk '{print $0"\t"$4"\t"$5}' 111 |sort -k7,8 |uniq -f 6 -u |cut -f1-6 >polploid_decrease.txt
## 多倍化后趋同
 sort -k5,6 22 |uniq -f4 -u >polploid_inrease.txt
## 多倍化后仍旧保守的
mv 333  polploid_unchange.txt

分析四倍体类保守转录本的组成

At与Dt是保守的，但是只在A2中保守，或者只在D5中保守；趋同进化

| At 与Dt 保守的6331 对isform中，有3166在二倍体内同样保守，1145对只在A2中出现，840对只在D5中出现，剩下1180对，其他

分析A2与D5中保守，但只在At与Dt一个亚基因组中保守的，并行分化

| A2与D5中保守的11823条转录本，在四倍体中有3166条同样保守，有2262条在At中保守，2380条在Dt中保守，4,015是其他类型

多倍化过程中两股力量同时推动着棉花的进化

虽然At与Dt比较中存在较少比例的保守isofrom，但是我们分别对A2与D5中保守的isofrom，与两个亚基因组比较后发现39%，4600多条只在At或者只在Dt中保守的isofrom；与此同时通过对At与Dt中已经鉴定到的保守的isform，我们发现有将近31%的isofrom只在A2或者只在D5中保守，更多的D基因组向A基因组方向进化

二倍体A2、D5合并后的isform数目与At、Dt合并后是否存在差异

  • At相比于A2 isform数减少

  • Dt相比于D5 isform数减少

  • At与Dt合并后相比于A2、D5合并后isform数减少

  • A2 与D5 间存在差异而在At、Dt中没有差异

保守的isform表达水平与非保守isform表达水平比，

A2 conserve vs A2 uniq | At conserve vs At uniq

保守的isoform表达水平更高

##获取每种比较中isform和对应的fpkm值
python ~/scripte/Alternative/module/FEST3/isformlevel/isformFPKM.py  -c conserve_isform.txt  -all isform_count.txt  -fpkm1 ../../geneExpress/stringtie/A2/t_data.ctab  -fpkm2 ../../geneExpress/stringtie/D5/t_data.ctab   -o A2_D5_conserveFPKM.txt  -uniq A2_D5_uniqFPKM.txt

四组同源基因中保守的isform表达水平的变化是否存在剂量效应

保守的isform与非保守的isform发生AS的频率

##获取发生剪切的isform
awk '{print $($8)}' ~/work/Alternative/result/Gr_result/CO41_42_result/11_AS/end_splice.txt3 |sort |uniq >splice_isform.txt
##提取每种类型中剪切的数目
awk '$1~/Ghir/{print $2}' A2_At_uniq_FPKM.txt |xargs  -I {} grep -E "{}$" ../TM1/splice_isform.txt |wc -l

比较

保守isform

发生剪切isform

不保守isform

发生剪切isform

A2 vs At

12868/12868

1445/1690

23680/13270

7232/4680

D5 vs Dt

12335/12335

990/1471

17868/14053

5305/4817

A2 vs D5

11823/11823

1121/841

21264/15536

6533/4750

At vs Dt

6331/6331

828/852

10286/10324

3187/3283

构建卡方检验表：AS与isform的保守性有无关系

保守的isform

不保守的isform

发生AS

没有发生AS

> table=matrix(c(11423,1445,16448,7232),nrow=2)
> chisq.test(table)

在多倍化后At、Dt亚组保守的isform数减少了，统计在两个二倍体多倍化后四倍体内亚组间的差异

  1. 在二倍体中保守的转录本

    1. 只在At中存在

    2. 只在Dt中存在

    3. At、Dt也存在

  2. 在二倍体中不保守的

    1. 只在At中存在

    2. 只在Dt中存在

    3. At、Dt都存在

在category目录中放了这几类

 python ~/scripte/Alternative/module/FEST3/fourHomologIsform.py  -c1 ../A2_vs_D5/conserve_isform.txt  -c2 ../A2_vs_At/conserve_isform.txt  -c3 ../At_vs_Dt/conserve_isform.txt  -c4 ../D5_vs_Dt/conserve_isform.txt  -u1 ../A2_vs_D5/A2_uniq.txt  -u2 ../A2_vs_D5/D5_uniq.txt -o1 1 -o2 2 -o3 3 -o4 4 -o5 5  -A2 A2_uniq -D5 D5_uniq

比较保守的isform之间表达量是否存在差异

##脚本内对两个基因组总的表达水平进行了矫正
python ~/scripte/Alternative/module/FEST3/conserveIsformFPKM.py -c conserve_isform.txt  -fpkm1 ../../geneExpress/stringtie/A2/t_data.ctab  -fpkm2 ../../geneExpress/stringtie/TM1/t_data.ctab  -o hhhh
##获取差异倍数变化
awk -v fild1=5 -v fild2=6 -f ~/github/zpliuCode/script/log2.awk hhhh >hhhh2

isform水平的差异表达

比较

总数

上调表达

下调表达

A2 D5

11823

2341

2650

At Dt

6331

1114

1193

A2 At

12868

3273

2300

D5 Dt

12335

3182

2120

isform表达水平的剂量效应

保守的isform，与不保守的isform平均表达水平差异

特异性的isform

python ~/scripte/Alternative/module/FEST3/uniqisform.py -c conserve_isform.txt  -a isform_count.txt  -u1 A2_uniq.txt  -u2 D5_uniq.txt

参考

  1. wu-blast

  2. ab-blast

Previous重新下载原始数据进行比对Next01

Last updated