龙岩网站建设要多久,网站集群建设方案,wordpress 登录 手机版,jsp网站开发标准读取 reads#xff08;二者含义相同#xff0c;下文不做区分#xff09;1. ChIPseq reads 比对 在评估读取质量和我们应用的任何读取过滤之后#xff0c;我们将希望将我们的读取与基因组对齐#xff0c;以便识别任何基因组位置显示比对读取高于背景的富集。 由于 ChIPseq…读取 reads二者含义相同下文不做区分 1. ChIPseq reads 比对 在评估读取质量和我们应用的任何读取过滤之后我们将希望将我们的读取与基因组对齐以便识别任何基因组位置显示比对读取高于背景的富集。 由于 ChIPseq 读数将与我们的参考基因组连续比对我们可以使用我们在之前中看到的基因组比对器。生成的 BAM 文件将包含用于进一步分析的对齐序列读取。 2. 参考基因组生成 首先我们需要以 FASTA 格式检索感兴趣的基因组的序列信息。我们可以使用 BSgenome 库来检索完整的序列信息。对于小鼠 mm10 基因组我们加载包 BSgenome.Mmusculus.UCSC.mm10。 library(BSgenome.Mmusculus.UCSC.mm10)BSgenome.Mmusculus.UCSC.mm10 BSgenome.Mmusculus.UCSC.mm10 我们将仅使用主要染色体进行分析因此我们可能会排除随机和未放置的重叠群。在这里我们循环遍历主要染色体并根据检索到的序列创建一个 DNAStringSet 对象。 mainChromosomes - paste0(chr, c(1:19, X, Y, M))mainChrSeq - lapply(mainChromosomes, function(x) BSgenome.Mmusculus.UCSC.mm10[[x]])names(mainChrSeq) - mainChromosomesmainChrSeqSet - DNAStringSet(mainChrSeq)mainChrSeqSet mainChrSeqSet 现在我们有了一个 DNAStringSet 对象我们可以使用 writeXStringSet 来创建我们的 FASTA 序列文件来比对。 writeXStringSet(mainChrSeqSet, BSgenome.Mmusculus.UCSC.mm10.mainChrs.fa) 3. 索引创建 我们将使用 subread 背后的 subjunc 算法进行对齐。因此我们可以使用 Rsubread 包。在我们尝试比对我们的 FASTQ 文件之前我们需要首先使用 buildindex() 函数从我们的参考基因组构建一个索引。 buildindex() 函数仅采用我们所需的索引名称和要从中构建索引的 FASTA 文件的参数。 library(Rsubread)buildindex(mm10_mainchrs, BSgenome.Mmusculus.UCSC.mm10.mainChrs.fa, memory 8000, indexSplit TRUE) 请记住建立索引会占用大量内存默认情况下设置为 8GB。这对于您的笔记本电脑或台式机来说可能太大了。 4. 比对 4.1. Rsubread 我们可以使用 Rsubread 包将 FASTQ 格式的原始序列数据与 mm10 基因组序列的新 FASTA 文件进行比对。具体来说我们将使用 align 函数因为它利用了 subread 基因组比对算法。 myMapped - align(mm10_mainchrs, filtered_ENCFF001NQP.fastq.gz, output_format BAM, output_file Myc_Mel_1.bam, type dna, phredOffset 64, nthreads 4) 4.2. Rbowtie2 Bowtie 家族是最著名的对齐算法之一。我们可以使用 Rbowtie2 包访问 Bowtie2。QuasR 包允许访问原始的 Bowtie 对准器但它有点慢并且需要内存。 library(Rbowtie2) 与 Rsubread 一样Rbowtie2 包要求我们首先创建一个要对齐的索引。我们可以使用 bowtie2_build() 函数来完成此操作指定我们的 FASTA 文件和所需的索引名称。 bowtie2_build(references BSgenome.Mmusculus.UCSC.mm10.mainChrs.fa, bt2Index file.path(BSgenome.Mmusculus.UCSC.mm10.mainChrs)) 然后我们可以使用 bowtie2() 函数对齐我们的 FASTQ 数据指定我们新创建的索引、SAM 输出的所需名称和未压缩的 FASTQ。 我们需要先解压缩我们的 FASTQ。这里我们使用 remove is FALSE 设置来保持原始压缩的 FASTQ。 library(R.utils)gunzip(filtered_ENCFF001NQP.fastq.gz, remove FALSE)bowtie2(bt2Index BSgenome.Mmusculus.UCSC.mm10.mainChrs, samOutput ENCFF001NQP.sam, seq1 filtered_ENCFF001NQP.fastq) 由于 Rbowtie2 还输出 SAM 文件我们需要将其转换为 BAM 文件。我们可以使用 RSamtools 的 asBam() 函数来做到这一点。 bowtieBam - asBam(ENCFF001NQP.sam) 使用 Rbowtie2 时的一个重要考虑因素是其未压缩文件的输入和输出。在命令行上我们可以将输入流式传输到 Rbowtie2但在 R 中这不是一个选项。我们需要确保删除任何创建的临时文件SAM 和/或未压缩的 FASTQ以避免填满我们的硬盘。我们可以使用 unlink() 函数删除 R 中的文件。 unlink(ENCFF001NQP.sam) 4.3. 排序 和以前一样我们分别使用 Rsamtools 包 sortBam() 和 indexBam() 函数对文件进行排序和索引。生成的排序和索引 BAM 文件现在可以用于外部程序例如 IGV也可以用于 R 中的进一步下游分析。 library(Rsamtools)sortBam(Myc_Mel_1.bam, SR_Myc_Mel_rep1)indexBam(SR_Myc_Mel_rep1.bam) 广告 本文由 mdnice 多平台发布
