目录
一、比对
首先在ensmbl上下拟南芥基因组序列(fasta格式)文件和基因组注释文件(GTF格式)。一定都要在ensmbl上下载,NCBI上下载GTF是GFF2.2qualimap无法使用。
- 访问Ensembl数据库:打开你的浏览器,访问Ensembl数据库的植物部分,网址为:Ensembl Plants。
- 查找拟南芥:在首页上,你会看到一些常用的物种列表,包括拟南芥。如果拟南芥不在首页上,你可以点击“View full list of all Ensembl Plants species”来获取所有物种的列表,然后从中找到拟南芥。
- 进入拟南芥参考基因组页面:点击拟南芥的链接,进入拟南芥参考基因组的介绍页面。
- 下载FA文件:在拟南芥参考基因组页面上,找到下载链接或相关选项。通常,你可以找到“Download DNA sequence (FASTA)”这样的选项。点击该选项后,你会看到多个文件可供选择。为了获取完整的参考基因组文件,通常建议下载名为
*toplevel.fa.gz
的文件。这是参考基因组的完整FASTA格式文件。 - Index of /pub/plants/release-58/gtf/arabidopsis_thaliana (ebi.ac.uk)
https://ftp.ensemblgenomes.ebi.ac.uk/pub/plants/release-58/gtf/arabidopsis_thaliana/网站如上,命名有规律,点击gff3然后将网址gff3改为gtf即可找到gtf下载地址。
- 保存文件:下载完成后,将文件保存到你希望存储的位置,并解压(如果需要)以供后续使用。
开启共享文件夹功能,在windows系统中将下载的文件放入共享文件夹中,共享文件位置在/mnt/hgfs/fold name,fold name即在windows系统中共享文件夹的名称。
1.建立基因组索引
安装好 Bowtie2 后,打开终端,切换到包含你希望保存索引文件的目录,然后运行下述命令。
gunzip -c Arabidopsis_thaliana.TAIR10.dna.toplevel.fa.gz > tair10.fa
bowtie2-build /home/lumino/TEST1/TEST1.1/tair10.fa tair10
1.gunzip -c:解压 Arabidopsis_thaliana.TAIR10.dna.toplevel.fa.gz 文件,并将解压后的内容输出到标准输出(而不是直接写入到一个新文件或覆盖原文件)。
2.>:重定向操作符,它将左边命令(在这里是 gunzip -c)的标准输出写入到右边的文件(在这里是 tair10.fa)。
3.bowtie2-build 是用来建立索引的命令,/home/lumino/TEST1/TEST1.1/tair10.fa 是你的基因组 FASTA 文件的路径,tair10 是你希望为索引文件使用的前缀。执行这个命令后,会在当前文件夹下生成一系列以 tair10 为前缀的索引文件,如 tair10.1.bt2、tair10.2.bt2、tair10.3.bt2、tair10.4.bt2 以及 tair10.rev.1.bt2 和 tair10.rev.2.bt2。
2.作图
#确保bowtiealign子文件夹存在。如果不存在,你需要先创建它
mkdir -p /home/lumino/TEST1/TEST1.1/bowtie/bowtiealign
# 确保你已经使用bowtie2-build为tair10.fa创建了基因组索引
# 确保你在包含SRR3418005-filtered.fastq文件的目录中
# 使用bowtie2进行比对,并指定索引文件路径/bowtie和输出文件路径/bowtie/bowtiealign
bowtie2 -x bowtie/tair10 -U SRR3418005-filtered.fastq -S /home/lumino/TEST1/TEST1.1/bowtie/bowtiealign/SRR3418005.sam
-
-x tair10告诉Bowtie2使用名为tair10的基因组索引进行比对。-U SRR3418005-filtered.fastq指定了输入的单端FASTQ文件。 -S /home/lumino/TEST1/TEST1.1/bowtie/bowtiealign/SRR3418005.sam
则指定了输出SAM格式文件的路径和名称。
-
- 5069801 (20.61%) aligned 0 times:这意味着20.61%的reads(即5,069,801个)没有与参考基因组或任何其他序列对齐。这可能是因为它们包含太多的错误、它们来自未知的区域或它们是非编码的。
- 16907641 (68.72%) aligned exactly 1 time:68.72%的reads(即16,907,641个)只与参考基因组或其他序列对齐了一次。这是预期的主要对齐方式,因为它表示这些reads来自基因组中的独特位置。
- 2626679 (10.68%) aligned >1 times:10.68%的reads(即2,626,679个)与参考基因组或其他序列对齐了多次。这通常表示这些reads来自基因组中的重复区域或高度相似的区域。
- 79.39% overall alignment rate(较为重要的指标):这是所有与参考基因组或其他序列至少对齐一次的reads的百分比。在这里,79.39%的reads至少对齐了一次,这表示大多数reads都能与参考序列匹配。
3.SAM文件处理
以BAM格式存储比对节省空间,并且许多下游工具使用BAM格式,而不是SAM。
使用samtools将SAM文件转换为BAM。去报提前下载samtools。
#进入SAM文件bowtiealign文件夹
cd ./bowtiealign
#使用samtools view命令,并指定输出格式为BAM
samtools view -b -o /home/lumino/TEST1/TEST1.1/bowtie/bowtiealign/SRR3418005.bam /home/lumino/TEST1/TEST1.1/bowtie/bowtiealign/SRR3418005.sam
#排序BAM文件
samtools sort /home/lumino/TEST1/TEST1.1/bowtie/bowtiealign/SRR3418005.bam -o /home/lumino/TEST1/TEST1.1/bowtie/bowtiealign/SRR3418005_sorted.bam
#为排序后的BAM文件创建索引,创建一个bai文件,bam和bai常用于后续基因组分析
samtools index /home/lumino/TEST1/TEST1.1/bowtie/bowtiealign/SRR3418005_sorted.bam
4.比对质量评估
#确保已经安装了QualiMap软件
# 假设当前工作目录包含tair10.gtf文件
# 运行QualiMap来评估比对质量
qualimap bamqc -bam /home/lumino/TEST1/TEST1.1/bowtie/bowtiealign/SRR3418005_sorted.bam -outdir /home/lumino/TEST1/TEST1.1/qualimap_analysis -gff ./tair10.gtf -nt 4
-
Mapped reads: 79.39% 的 reads 被映射到了参考序列上,这是一个相对较高的映射率。这意味着大部分测序数据都能够与参考基因组对齐,从而可以用于后续的分析。
-
Unmapped reads: 20.61% 的 reads 没有被映射,这通常是由于这些 reads 含有过多的测序错误、来自基因组的未知区域或属于非编码区等原因。这个比例并不是特别高,但也可能意味着一部分信息丢失。
-
Mapped paired reads: 配对 reads 的映射数量为 0,这意味着您提供的测序数据是未配对的,或者是单端测序数据。如果是预期为配对测序的实验,这可能是数据准备或处理过程中的一个问题。但如果实验设计就是单端测序,那么这个数据就是正常的。
-
Secondary alignments: 没有次级对齐的 reads,这通常表示数据中的 reads 都是唯一的,没有来自重复区域或高度相似区域的 reads。
-
Read length: 所有 reads 的长度都是 89 bp,这是一个相对标准的测序长度,适用于多种应用。长度的一致性也表明测序过程中没有显著的长度偏差。
-
Clipped reads: 没有 reads 被裁剪,这意味着在数据预处理阶段没有显著的序列质量问题需要修剪。
-
基于上述分析,我们可以认为这些测序数据在映射率和读取长度方面表现良好,但有一定比例的 reads 未映射,可能需要注意这部分信息的缺失对后续分析的影响。然而,数据质量的全面评估还需要结合更多的质量指标,如碱基质量分数、GC含量分布等,这些可以通过专业的质量控制工具(如FastQC)来检查。
综上所述,从提供的数据来看,这些数据的质量在映射率和读取长度方面看起来是合理的,但还需要进一步的质量控制分析来确认其完整的质量状况。