AS保守程度的统计

针对各个基因组之间比较得到的保守性，将基因对分类

甲基化程度的检测

将多重校验的结果做成bed文件

awk '{print $1,$2,$2,$6}' OFS="\t" fdr.file >fdr.bed

制作基因组的bed文件

Ghir_A11    48563862    48571977    Ghir_A11G021610
Ghir_A11    48086239    48087160    Ghir_A11G021550
Ghir_A11    48315032    48323201    Ghir_A11G021570

使用intersectBed取交集

• -a 基因坐标文件

• -b 甲基化bed文件

• -loj 进行左连接，如果没有检测到胞嘧啶则会为NULL

~/software/bedtools2-2.29.0/bin/intersectBed -loj -a ../TM-1_gene.bed -b CpG_context_TM1_qs_rep1_count.txt_binom.txt_fdr.bed >3

计算甲基化程度

由于测序深度会影响一个基因甲基化程度的检测，所以甲基化程度的计算方法

被检测到的甲基化的C碱基 / 被测到的所有C碱基

awk '$8=="."{print $4"\t0"}$8!="."&&$8<=0.05{a[$4][0]+=1}$8!="."&&$8>0.05{a[$4][1]+=1}END{for(i in a){print i"\t"a[i][0]/(a[i][0]+a[i][1])}}' 1 >2

检测重复性

计算pearson相关系数

data1 <- read.table("CHG_context_TM1_qs_rep1_count.txt_binom.txt_fdr.bed_end_sorted")
data2 <- read.table("CHG_context_TM1_qs_rep2_count.txt_binom.txt_fdr.bed_end_sorted")
cor.test(data1$V2,data2$V2,method = c("pearson"))

画图

重复情况检测

参考

相关性 https://www.jianshu.com/p/b76f09aacd9c