DNA methylation plays a crucial role in most organisms. Bisulfite sequencing (BS-Seq) is an approach that can be used to obtain quantitative cytosine methylation levels at a genome-wide scale and single-base resolution. Besides, parental alleles in haploids might exhibit different methylation patterns, which can lead to different phenotypes and even different therapeutic and drug responses to diseases. A comprehensive database of high-throughput BS-Seq data and allele-specific DNA methylation (ASM) results from multiple species was not available prior to this study.

Here, we constructed The Allele-Specific DNA Methylation Databases (ASMdb), aiming to provide a comprehensive resource and a web tool for showing the DNA methylation level and differential DNA methylation in diverse organisms, including 47 species, such as Homo sapiens, Mus musculus, Arabidopsis thaliana, and Oryza Sativa, with 5998 GEO samples (4400 BS-Seq data and 1598 RNA-Seq data).

1. The trimming of low-quality reads and artificial sequences was performed with fastp.
2. The clean reads were mapped to reference genomes using BatMeth2-align, and the SAM files were converted to the BAM format with SAMtools.
3. DNA methylation calling was performed with calmeth. After sequences with a map quality score lower than 20 were filtered out, the cytosine sites with a coverage greater than 5 were considered effective methylation sites.

1. The alignment file for each sample was obtained using the BatMeth2 software package, and the alignments were then sorted using SAMtools.
2. SNP calling was performed with BSsnpcall with the default parameters.

1. Calculation of the DNA methylation level and DNA methylation alignment were performed using the BatMeth2 software package.
2. ASM was identified by MethHaplo with the default parameters.

1. The raw reads from the RNA-Seq data were first trimmed as paired-end reads using fastp with the default parameters to remove the adaptors and low-quality reads.
2. The clean reads were mapped to the reference genomes using Hisat2, and SAMtools was then used to sort the BAM file.
3. The SNP information was derived from the BS-Seq analysis results corresponding to the RNA-Seq data. ASEGs were detected by ASEQ..

1. DNA methylation BED files containing the coverage and methylation levels were analyzed using BatMeth2.
2. Each DNA methylation sample was divided into partially methylated domains (PMDs), lowly methylated regions (LMRs) and unmethylated regions (UMRs) according to the methylation level using MethylSeekR.
3. Furthermore, we removed the 'N' gap region in the genome, and the remaining region is called highly methylated domains (HMDs)
4. The details of the scripts are shown in below:

   #MethylSeekR
   install.packages("BiocManager")
   BiocManager::install("MethylSeekR")
   library(MethylSeekR)
   library(BSgenome)
   library(parallel)
   library(rtracklayer)
   
   BiocManager::install("BSgenome")
   library(BSgenome)
   #available.genomes()
   BiocManager::install("BSgenome.Osativa.MSU.MSU7")
   library("BSgenome.Osativa.MSU.MSU7")
   build <- get("BSgenome.Osativa.MSU.MSU7")
   #also can bulid by ourself..
   #library(BSgenome.Astyanax.mexicanus)
   #build <- get("BSgenome.Astyanax.mexicanus")
   Args <- commandArgs()
   faifile<-Args[6] ##faifile
   pref<-Args[7] ## pref
   cpgisland<-Args[8] ##cpgislad
   trainchr<-Args[9] ##train chr
   fai <- read.delim(faifile, header=F)
   chromosome_lengths <- fai$V2
   names(chromosome_lengths) <- fai$V1
   d <- readMethylome(paste(pref, "_loci.CG.txt.methylseekr", sep=""), chromosome_lengths)
   snps.gr <- readSNPTable(paste(pref, ".snp.methylseekr", sep="") , chromosome_lengths)
   ## removr snp
   meth.gr <- removeSNPs(d, snps.gr)
   sample = pref
   num.cores=4
   pdf(paste0(sample, ".pdf"), width=14)
   PMDsegments.gr <- segmentPMDs(m=meth.gr, chr.sel=trainchr, seqLengths=chromosome_lengths, num.cores=num.cores)
   ##
   CpGisland <- import(cpgisland, format="bed")
   write.csv(calculateFDRs(meth.gr, CpGisland, PMDs=PMDsegments.gr, num.cores=num.cores)$FDRs, paste0(sample,".PMD.stats"))
   write.csv(calculateFDRs(meth.gr, CpGisland, num.cores=num.cores)$FDRs, paste0(sample,".stats"))
   ### UMR LMR
   PMD.UMRLMR <- segmentUMRsLMRs(m=meth.gr, meth.cutoff=0.5, nCpG.cutoff=5, num.cores=num.cores, myGenomeSeq=build, seqLengths=chromosome_lengths, PMDs=PMDsegments.gr)  
   # Repeat without filtering PMDs
   UMRLMR <- segmentUMRsLMRs(m=meth.gr, meth.cutoff=0.5, nCpG.cutoff=5, num.cores=num.cores, myGenomeSeq=build, seqLengths=chromosome_lengths)
   dev.off()
   # Export to bed/bigbed
   tmp <- list("PMD"=PMDsegments.gr[PMDsegments.gr$type=="PMD"],
               "UMR"=UMRLMR[UMRLMR$type=="UMR"],
               "LMR"=UMRLMR[UMRLMR$type=="LMR"],
           "PMD.UMR"=PMD.UMRLMR[PMD.UMRLMR$type=="UMR"],
           "PMD.LMR"=PMD.UMRLMR[PMD.UMRLMR$type=="LMR"])
   chrom.sizes <- data.frame(seqlevels(build), seqlengths(build))
   write.table(chrom.sizes, paste0(sample, ".chrom.sizes"), quote=FALSE, sep="\t", row.names=FALSE, col.names=FALSE)
   for (j in names(tmp)) {
       bed.file.name <- paste(sample, j, "bed", sep=".")
       export(tmp[[j]], bed.file.name, format="bed")
       system(paste0("cut -f1,2,3 ", bed.file.name, " > tmp; /public/home/qwzhou/software/EMBOSS-6.6.0/bin/bedToBigBed tmp ", sample, ".chrom.sizes ", gsub(".bed$", ".bb", bed.file.name), "; rm tmp"))
   }
   

• The methylome browser is built on JBrowse (http://jbrowse.org; a fast, embeddable genome browser built completely with JavaScript and HTML5), click here for information on how to use it.