> For the complete documentation index, see [llms.txt](https://zpliu.gitbook.io/booknote/llms.txt). Markdown versions of documentation pages are available by appending `.md` to page URLs; this page is available as [Markdown](https://zpliu.gitbook.io/booknote/ke-bian-jian-qie/di-san-ci-dui-as-jin-hang-tong-ji/02.md).

# 02

## 保守的 isforms鉴定

使用ab-Blast进行鉴定

#### 使用xdformat建库

* `-n`指定序列类型
* -o 输出数据库前称
* 最后跟fasta文件

```bash
xdformat  -n  -o test test.fa
```

#### 进行比对

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

* 库前称
* 序列文件
* `-o`输出文件

```bash
~/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的情况

```bash
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的长度都接近

```bash
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的

```bash
##转换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的总数的平均值

```bash
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](https://s1.ax1x.com/2020/06/23/NN2lIP.png)

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

```bash
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多倍化后，在两个亚基因组间共享

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

```bash
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直接比对的时候没有鉴定出保守性，但是通过两个基因组间的关系间接证明了保守性，这个先去掉，原始文件先存着；到时候设计引物的时候在看看

```bash
## 多倍化后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](https://s1.ax1x.com/2020/07/02/NqjSq1.png)

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

* At相比于A2 isform数减少
* Dt相比于D5 isform数减少
* At与Dt合并后相比于A2、D5合并后isform数减少
* A2 与D5 间存在差异而在At、Dt中没有差异

![NDvp59.png](https://s1.ax1x.com/2020/06/26/NDvp59.png)

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

> A2 conserve vs A2 uniq | At conserve vs At uniq
>
> 保守的isoform表达水平更高

```bash
##获取每种比较中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](https://s1.ax1x.com/2020/06/27/NyQLfH.png)

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

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

```bash
##获取发生剪切的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 |           |            |

```r
> 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`目录中放了这几类

```bash
 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之间表达量是否存在差异

```bash
##脚本内对两个基因组总的表达水平进行了矫正
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](https://s1.ax1x.com/2020/06/20/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

```bash
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](https://www.advbiocomp.com/blast/obsolete/)&#x20;
2. [ab-blast](https://blast.advbiocomp.com/licensing/personal.html)