设计引物

之前看用RNA-seq去看转录本的表达量的的时候，由于模板浓度太低了，跑不出发生AS的那条带；因此使用PacBio的read数目来进一步筛选发生剪接的转录本。

#A2 vs At
evolution4/A2_vs_AT/Isoform_DomainFPKMChange/IsoformDomain_readCount.txt
##查看read数目和基因数目
grep -v "Nochange"  IsoformDomain_readCount.txt |cut -f1,3,8,9|head
#D5 vs Dt

使用发生AS的基因和productive transcript发改变的基因取交集，同时满足PacBio reads数目大于10

grep Nochange -v conserve_specialIsoformDomain.txt |cut -f3|cat - AS_gene.txt |sort |uniq -d|xargs  -I {} grep {} ../../../geneFPKM_readsCount/TM1_geneFPKM_read.txt

比较不同棉种中productive转录本序列是否发生较大的缺失

当提取两个棉种中主要表达的转录本进行，预测它们的CDS的时候；有时候它们的cDNA仅仅只相差个几bp，都可能会影响后续CDS的预测；因此还要进行进一步筛选，去比较两个cDNA是否存在较大的indel

脚本用法: https://github.com/BiocottonHub/zpliuCode/tree/master/transcriptSV#transcriptsv

##鉴定祖先和后代中主要表达的转录本是否存在PAV
python ~/github/zpliuCode/transcriptSV/transcriptSvs.py -homolog ~/work/Alternative/result/homologo/homologGene/Result/A2_vs_At_collinearity.txt -fasta1 ~/work/Alternative/result/Ga_result/CO11_12_result/collapse/mappingTogene.fa  -fasta2 ../../collapse/mappingTogene.fa -RNAseq1 ../A2_PacBio.txt -RNAseq2 ../TM1_PacBio.txt -out ./A2_vs_At_productive_isoform.txt

是否存在AS导致的productive 转录本发生改变。

这里主要表达的转录本还是看full-read 数目是不是最多的。

1. 首先提取同源基因中主要表达的转录本

2. 比较两者是否存在缺失。

3. 比较两者CDS序列是否存在差异

存在以下几种情况：

1. 两个转录本CDS序列一致，而转录本存在结构变异；也就是AS影响的是OUTR区域。

2. 两个转录本CDS序列不一致，并且转录本存在结构变异；说明可能是AS影响了编码区

3. 两个转录本的CDS序列一致，同时转录本不存在结构变异；说明productive 转录本是保守的。

4. 两种转录本的CDS序列不一致，而不存在结构变异，这就可能由于MUMmer只能鉴定比较长的PAV

例如Gorai.010G219800与Ghir_D06G020630相比存在一个主要表达转录本的差异，在D5、At中第一个内含子都是组成型的剪接，而在Dt中第一个外显子持续性的保持外显子的状态；这还是个热激蛋白。

从这个例子中，我们不难发现有的时候D5、Dt可能只转录出了一种类型的转录本，因此没能检测到可变剪接的存在，但将D5与Dt转录出的转录本进行比较，我们可以发现它们的转录本类似的发生了AS。

TM-1

A2

D5