表观遗传在AS中的作用2

比较constitutive 事件与AS事件在CG甲基化程度上的差异

参考 DNA-methylation effect on cotranscriptional splicing is dependent on GC architecture of the exon–intron structure

统计某一个位置出现CG的次数，以及CG被甲基化的比例

##统一保留内含子第一个和最后一个坐标
grep RI ~/work/Alternative/result/Gh_result/CO31_32_result/evolution2/A2_AS.txt |grep PB|cut -f3,4 >1
python ~/github/zpliuCode/Isoseq3/04ASconserved/extractAScoordinate.py  1 1 2
##RI坐标
awk '{OFS="\t";print $1,$2+1,$3,$4,$6}' 2  >RI.bed
##SE坐标
awk '{OFS="\t";print $1,$2,$3+1,$4,$6}' 2  >SE.bed
##对于组成型的exon和intron，根据基因的链，赋予对于的链,对于超出原有注释范围的不管
~/software/bedtools2-2.29.0/bin/intersectBed -a constitutive_intron.bed -loj -b ../TM1_gene.bed |awk '$4!="."{OFS="\t";print $1,$2,$3,$7,$9}' >1 
mv 1 constitutive_intron.bed    
##分割文件

##提取每个坐标处的是否甲基化
python ~/github/zpliuCode/Isoseq3/07methylation/methylationByLocation.py  ~/work/Alternative/data/Ga_genome/G.arboreum.Chr.v1.0.fasta  /data/cotton/zhenpingliu/QingxinSong_GB_DNAmethlation/A2/Rep1/02deduplicate_methylation/CpG_fdr_sorted.bam 1 3
##合并所有的任务的结果
python  ~/github/zpliuCode/Isoseq3/07methylation/gatherAlleventLocal.py conExon_methylation.txt
##由于基数太大可能会造成差异，只对AS基因进行比较

比较AS基因和非AS基因的甲基化差异程度

方法参考： Asymmetric epigenome maps of subgenomes reveal imbalanced transcription and distinct evolutionary trends in Brassica napus

基因区域被分为了40个bin，以及基因上下游2kb分同样各分为40个bin
每个bin的甲基化程度用（甲基化的read/测到的总的read）进行衡量

##提取上游的2kb坐标
 ~/software/bedtools2-2.29.0/bin/flankBed -b 2000 -i TM1_gene.bed  -g TM1_chromsome.bed |awk '$NF=="+"&&NR%2==1{print $0}$NF=="-"&&NR%2==0{print $0}' >Upstream2K_gene.bed

 ##做成40个bin的窗口
 ##统一用远离基因的一端做起始bin，区分正负链
awk '$NF=="-"{print $0}' Upstream2K_gene.bed |~/software/bedtools2-2.29.0/bin/windowMaker -b - -n 40 -i srcwinnum -reverse >Upstream2k_40bin.bed
#正链则不需要反向
awk '$NF=="+"{print $0}' Upstream2K_gene.bed |~/software/bedtools2-2.29.0/bin/windowMaker -b - -n 40 -i srcwinnum >>Upstream2k_40bin.bed
##提取基因某一个区域的甲基化程度
python ~/github/zpliuCode/Isoseq3/07methylation/extractMethylationRegion.py  gene_40bin.bed /data/cotton/zhenpingliu/QingxinSong_GB_DNAmethlation/TM1/Rep1/02deduplicate_methylation/CpG_fdr_sorted.bam gene_40bin_CpG.txt 
##合并所有的bin，求平均的甲基化水平
cat Downstream/*|sed 's/_/\t/g' |awk '{a[$2][1]+=$3;a[$2][2]+=1;}END{for(i in a){print i "\t"a[i][1]/a[i][2]}}'|sort -k1,1n|less
##区分AS基因和非AS基因的甲基化程度

对AS与非AS基因进行显著性检验

genome

CpG

CHG

CHH

TM1

0.00159

0.5531

0.2222

1.303e-05

0.1628

0.6058

0.0001102

0.8664

0.1216

conExon与SE区域的平均甲基化程度没有差异，而RI与conIntron之间存在差异

分析不同基因组的事件间序列的保守程度

在筛选组成性内含子的时候，用了merge这一步骤导致很多intron被合并了，所以对于差异比较大的片段还是使用k-mer的坐标比较准确一些

##筛选保守的AS事件的坐标
python ~/github/zpliuCode/Isoseq3/07methylation/filterConstitutiveintronByKmer.py  ~/work/Alternative/result/Gh_result/06ASconserved/test2/A2_vs_At/AS_kmer.txt A2RI_2_Atintron.txt  A2RI_2_Atintron_filter.txt 

##统计序列变异情况
python ~/github/zpliuCode/Isoseq3/07methylation/sequence/sequenceAligment.py ~/work/Alternative/data/Ga_genome/G.arboreum.Chr.v1.0.fasta ~/work/Alternative/data//Ghirsutum_genome_HAU_v1.0/Ghirsutum_genome_HAU_v1.0.fasta A2_CI.txt At_RI.txt Donor_variant.txt Acceptor_variant.txt SequenceIdentity

比较序列完全保守的坐标之间CpG甲基化程度的差异

##序列完全保守的区域，甲基化程度的差异
python ~/github/zpliuCode/Isoseq3/07methylation/sequence/extractSequenceMthylation.py  /data/cotton/zhenpingliu/QingxinSong_GB_DNAmethlation/D5/Rep1/02deduplicate_methylation/CpG_fdr_sorted.bam /data/cotton/zhenpingliu/QingxinSong_GB_DNAmethlation/TM1/Rep1/02deduplicate_methylation/CpG_fdr_sorted.bam D5RI_2_DtCI_identity.txt  1

分析同源基因间保守的胞嘧啶坐标

根据同源基因的坐标筛选保守的C坐标
获取区域内总的胞嘧啶个数
~~同源基因中CG满足40%的覆盖度~~

##筛选胞嘧啶保守的位置
python ~/github/zpliuCode/Isoseq3/07methylation/sequence/conservedCytosine.py ~/work/Alternative/data/Ga_genome/G.arboreum.Chr.v1.0.fasta ~/work/Alternative/data/Ghirsutum_genome_HAU_v1.0/Ghirsutum_genome_HAU_v1.0.fasta ../coordinate/A2_homolog_gene.bed  ../coordinate/At_homolog_gene.bed  A2_At_conservedCystoin.txt
##获取区域内CG碱基对数目，以及满足read数大于3的CG碱基对数目
python ~/github/zpliuCode/Isoseq3/07methylation/sequence/geneWithMethylatedCytosinesCount.py  CpG1.bam CpG2.bam geneLocation.txt genome.fa
##筛选同源基因长度CG满足40%的覆盖度

对保守的C进行DmC分析

同源基因CG甲基化的差异程度：
DmCG数目/同源基因保守的CG碱基对数目

##对每个保守的胞嘧啶位点进行差异分析
python ~/github/zpliuCode/Isoseq3/07methylation/sequence/DmCG.py conservedCytosine.txt  conservedCytosineRatio.txt CpG_Bam 
CpG_Bam  CpG_Bam  CpG_Bam 
##获取每对同源基因CG甲基化的差异程度
awk '{a[$5][1]+=1}$11<=0.01{a[$5][2]+=1}END{for(i in a){print i"\t"a[i][1]"\t"a[i][2]}}' conservedCytosineCpG.txt|awk '$3==""{print $1"\t"$2"\t0\t0.0"}$3!=""{print $0"\t"$3/$2}'

比较

保守的胞嘧啶

保守的CG位点

CG差异甲基化

A2 vs At

10,613,071

1,100,398

173319 (15.8%)

D5 vs Dt

12,278,190

1,272,223

175363 (13.8%)

A2 vs D5

10,180,352

1,086,003

207524 (19.1%)

At vs Dt

11,622,245

1,226,762

210935 (17.2%)

基因的DmCG与基因AS差异程度的比较

基因AS的差异程度主要比较的是两个同源基因都存在AS的基因对

hyper conserved AS ： AS保守程度大于0.8

hypo conserved AS ： AS保守程度小于0.4

AS高度保守的基因相比于低度保守基因，DmCG/CG 比例更小
序列水平的变异

# result/Gh_result/06ASconserved/02conservedAS/quantifyASconserveRatio

类别

p-value

A2 At

6.92e-06

A2 D5

2.2e-16

D5 Dt

2.2e-16

At Dt

2.2e-16

基因AS数目的差异与DmCG/CG

AS数目差异越大，DmCG差异也越多

##获取同源基因AS数目的差值
awk '{print $1"\t"$2"\t"sqrt(($3-$4)*($3-$4))}' D5_vs_Dt_AScount_DmCG.txt

源基因序列发生变异附近的motif 富集分析

统一用供体开始的坐标计算变异位点
提取变异位点左右各10bp的序列进行motif富集分析

##提取变异位点左右各10bp的序列进行motif富集分析
python ~/github/zpliuCode/Isoseq3/07methylation/sequence/extractNearSequenceVarant.py genome1.fa genom2.fa sequenceVariant.txt motif20bp.fa
##进行motif 富集分析

A2RI : AAGMCCTGSKAAYMTG
A2CI : GGMCHAGCCTRTCACT
AtRI : TCTTYCYTTYHC
AtCI : AARARGAARRWRRAAR

D5RI : VCCHYCYCCYCCCCC
D5CI : NDCCCSMCCCCCCCC
DtRI : MTTTTTTTTTTT
DtCI : YGTTYGRGRMAWGAA

AS差异数目

p-value

0 vs 1

2.2e-16

1 vs 2

2.2e-16

2 vs 3

2.2e-16

3 vs 4

2.2e-16

4 vs 5

2.2e-16

举个例子

特异性的AS事件中，存在甲基化差异的事件

##获取亚组特异性的剪接事件

类型

总特异性事件数

存在甲基化差异事件数

A2_vs_At A2

15126

6789 (44.9%)

A2_vs At At

10458

4840 (46.3%)

A2_vs D5 A2

15339

7706 (50.2%)

A2_vs D5 D5

15361

7569 (49.3%)

D5_vs_Dt D5

14525

7321 (50.4%)

D5_vs_Dt Dt

9781

4981 (50.9%)

At_vs_Dt At

10803

5988 (55.4%)

At_vs_Dt Dt

10922

5893 (54.0%)

Total

102,315

51,087

evm.TU.Ga03G2363基因

Previous全长转录本数据的统计 Next03表观遗传与可变剪切

Last updated 3 years ago