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

NN2lIP.png

二倍体到四倍体，多倍化过程中，存在保守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基因组方向进化

NqjSq1.png

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

• At相比于A2 isform数减少

• Dt相比于D5 isform数减少

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

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

NDvp59.png

保守的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

NyQLfH.png

四组同源基因中保守的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

Nla5uR.png

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